OTS PRICE 



NsG-'or-Lf 

N 63 19946 

Ca^UjLs / 



College of Engineering 

RESEARCH DIVISION 

University Heights, New York 53, N. Y. 


GEOMAGNETIC IMPLICATIONS OF THE 
COSMIC RAY DIURNAL VARIATION, 
AND ITS SPECTRUM 

by 

Donald E. Cotten 



Final Project Report 


Prepared for 

U.S. Army Signal Research and Development Laboratory 

and 

National Aeronautics and Space Administration 

March 1, 1963 


NEW YORK UNIVERSITY 
College of Engineering 
Research Division 


GEOMAGNETIC IMPLICATIONS OP THE 
COSMIC RAY DIURNAL VARIATION 
AND ITS SPECTRUM 
by 

Donald E. Cotten 


Pinal Project Report 


Approved 


Arthur Beiser 
Project Director 


The work described herein was sponsored by the A J^anced Research 
Projects Agency under Order Nr. 163-61, through the U.S Army 
Signal Research and Development Laboratory under Contract Nr. 

DA 36-039-SC-87171, and by the National Aeronautics and Space 
Administration under Grant No. NsG 108-61. 


c.? 


Acknowledgements 


The author wishes to thank Dr. Arthur Beiser for 
suggesting the problem, and to express his gratitude to 
Drs. Arthur Beiser, Serge A. Korff, Wallace Arthur, Walter 
Kertz, Wolfgang Ramm, Jurgen Untiedt and David Stern, and 
also Thomas Kelsall, for many stimulating and helpful dis- 
cussions. Thankb are also due to Emmanuel Mehr and Abe 
Lucas for their assistance in computer programming, to 
Robert H. Kress for neat tracing of some of the figures 
in this report, to Joan Santora for her conscientiousness 
in typing and to my wife Evelyn Cotten for her assistance 
in preparation of this manuscript and her patience during 
the course of the research. 


ii 


ABSTRACT 


//? * 4 


When the amplitudes n 1 of the diurnal variation (CRDV) in 
counting rates of all northern and equatorial neutron monitors 
are averaged over groups of months within the IGY, and some pairs 
of stations with close values of Quenby-Webber cutoff rigidity P c 
are averaged together to smooth the data, then two peaks of n 1 (P c ) 
persist at approximately P c =l and 4 BV/C. An expression for n^ 
wherein the impact zones as calculated by Kelsall for various 
rigidities P are expressed as geomagnetic colatitude-dependent 
integration limits P Q , ...P^, 4s used to solve for the CRDV differ- 
ential spectrum (4/3)k(P)cos 0 at those limits without any prior 


assumption as to the form of k(P). Here k(P) is a fraction of the 
isotropic cosmic ray flux at infinity and 0 is the assymptotic 
longitude with respect to an axis of anisotropy. If the equatorial 
data is smoothed, then k(P) = 0.0039 ± 0.002 satisfies the equation 
at 6 BV/C, and the double peak of n 1 (P (J ) leads to a single peak 
of k(P) at 3.8 <^<[5.7 BV/C, with maximum k(P) = 0.01 to 0.06 at 
P = 5.4 BV/C. At P<3.8 BV/C, k(P) = 0. Contours are presented 
showing the CRDV amplitude and phase, and the amplitude of the 
semi-diurnal variation, as functions of latitude and longitude. 
These maps display relative extremums which correspond to extre- 
mums of geomagnetic field maps, and which indicate that the geomag- 
netic dipole and quadrupole moments as measured at the earth's 


surface significantly affect the anisotropic part of the cosmic 


radiation. An explanation is given for the fact that the semi- 
diurnal variation has appreciable amplitude only along part of 
the magnetic equator. The cosmic ray anisotropy introduced by 

the solar system is found to come from a direction 75° to 100 
east of the sun. Small CRDV components dependent on Greenwich 

time and sidereal time are discussed. 
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I. Introduction 

A. 

A description of the magnetic field at large distances 
above the earth's surface could be obtained completely from an 
infinite number of magnetometers distributed over the entire surface 
if not for the existence of unknown electric currents near the earth. 
This failure makes it desirable to measure the geomagnetic field far 
above the surface. Nature provides us at all times with a vast 
number of charged cosmic ray particles which have probed a large 
portion of the outer field. Their orbits and points of impact 
upon the atmosphere near the earth's surface depend on the inter- 
planetary magnetic field, the distant geomagnetic field, and partly 
on the nearby field. Collisions with atoms or molecules of our 
atmosphere produce secondary protons, neutrons, mesons, electrons, 
and photons, any of which may be detected at the earth's surface. 

The counting rates of the detectdrs can be related to the primary 
cosmic ray flux . 1 

In order to gain information about magnetic fields 
these counting rates must be compared over a large region of the 
earth, preferably the entire earth. Such comparison can only be 
made for comparable percentage changes in counting rate which 
occur over large areas. Forbush-type decreases (FD) of a few per- 
cent in the cosmic radiation occur at irregular intervals a few 
days or weeks apart and provide one means of such comparison. 

Solar flares produce large increases which permit comparison at 
even less frequent intervals. These events provide occasional 
information regarding the kind of disturbance of the interplanetary 
or geomagnetic field which occurred during the event. In order 

1 Quenby, J. J. and Webber, W. R., Phil . Mag. 4, 657 (1959). 


to monitor the undisturbed magnetic fields by using cosmic rays it 
is necessary to study a change In their counting rate which is world- 
wide and which occurs regularly on geomagnetically undisturbed days. 
The diurnal variation In cosmic ray nucleon counting rate (CRDV) 
permits such comparison between stations, provided that a suitable 
theory is used to relate the observed CRDV to the magnetic fields 
or other geophysical phenomena. This paper presents a new theory 
which accounts for all of the general features of the CRDV and also 
several of its heretofore unnoticed dnd unexplained details. 


I. B. Interplanetary and Geomagnetic Fields 

The geomagnetic field is specified at all points on the 
earth’s surface and in nearby space by the spherical harmonic co- 
efficients^* ^ for a magnetic scalar potential chosen to fit magne- 
tometer measurements at a limited number of fixed or mobile observing 
stations on the surface or in airplanes. 2 Satellite-borne magne- 
tometers have recently carried the measurements out to several earth 

4 5 

radii within limited regions of latitude, longitude and altitude. * 
Scalar potential analysis is only valid in the region which does not 
include currents such as in the Van Allen particle belts or an 
equatorial ring current beyond that.^ In that analysis the dipole 
terms dominate to give a dipole moment of 8.1 x 10 gauss cm 
oriented at 12° with the rotation axis. The quadrupole moment 


^Vestine, E. H., Transactions A.G.U. 4l , 4 (i960). 

3Finch, H. F. and Leaton, B. R., Roy. Astron.Soc., Mon. Not., 
Geophys. Suppl. £, 314 (1957)* 

^Cain, J. C., Shapiro, I. R., Stolarik, J. D. Heppner, J. P., 

NASA Report X-611-62-128, Aug. 21, 1962, Jour. Geophys. Res. 67, 

5Heppner^ 9 j?^P., Ness, N. F., Scearce, C. S. and Ski liman, T. L. 

NASA Report X-6 62-125, Jour. Geophys. Res. 68, 1> (1963). 

6 Akasofu, S-I, and Chapman, S., Jourv. Geophys. Res. 66, 1321-1850, 
(1961). 
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gives a field oriented sc that its maximum strength is attained at 
four "poles" situated near the dipole equator. Its average strength 
at the surface is 1 % of that of the dipole field, ^ and above the surface 
it falls off as r- 4 while the dipole field falls off as r' 3 , where r 
is the distance from the earth's center. In the theory of many geo- 
physical phenomena the earth's field is regarded as consisting only 
of the field of a magnetic dipole located at the center of the earth. 
When a slightly more detailed model of the field is required a dipole 
of the same strength is located eccentrically so as to best fit the 
strength of the dipole plus (rotation) axial quadrupole terms. 8 More 
detailed models for the field include the remaining four quadrupole 
terms and then multipole terms of higher order. Calculations for 
some phenomena, as for example, the guiding of hydromagnetlc waves 
bythe actual field lines, or the location of conjugate points at the 
two ends of a field line, require the use of as many multipole co- 
efficients as are available. The earth ring currents also make con- 

Q 

tributions to the observed geomagnetic field. 

Some measurements have been made of the strengths of inter- 
planetary magnetic (IP) fields, 5 and various models of the IP field 
have been proposed to explain the observed strength and infer its 
direction. One of the more recent proposals is that the plasma 
called the solar wind continuously expanding from the solar corona 
carries with it magnetic field lines bent in an Archimedian spiral 
and co-rotating with the sun. The direction of the IP field lines 
in the vicinity of the earth is not yet established. 

W. R., Phil. Mag. 4, 90, (1959). 

Atmos. Elec. 4l, 225, (1936) 

., and Cain, J. C., J. Geophys. Res. 

10 Dessler, A. J., Ahluwalia, H. S. and Gottlieb, B., J. Geophys. 

Res. 61, 3553, (1962). 


7Quenby, J. J. and Webber, 
OBartels, J. Terr. Mag. and 
9Chapman, S., Akasafu, S. I 
66, 4013, (1961). 
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Interaction between the solar wind plasma and the geo- 
magnetic field produces a geomagnetic envelope 11 ' 12 * ^ich that the 
geomagnetic field is compressed somewhat and Increased in strength 
on the solar-windward side or head, and extended and decreased in 
strength on the opposite side, or tall, of the geomagnetic field. 


I. C. Previous Studies of the Cosmic Ray Diurnal Variation 

Early in the study of cosmic radiation (CR), indications 

were found of a very small CRDV with period of one sidereal day 121 

22 

and of one mean solar day 15-21 , and a CRDV theory was put forth. 

23 24 

Altitude dependence of the observed CRDV was small. ' Attempts 

were made 21 ' 25 to attribute the CRDV to a CR anisotropy produced 
by a solar magnetic moment of 10 8 ^ gauss cm 8 , or to a beam of 
particles of CR energy coming from the sun. 26 ' 2 ? Deflections by 
the geomagnetic dipole model field were calculated 25 ' 28-31 f or 


H Hurley, J., Doctoral Thesis, New York University, (1961). 
12 Hurlev, J., "Interaction Between the Solar Wind and the Geo- 

' 9 •• i / — « -m 1 T 1 \ 


J-^Compton, A. H. and Getting, I. A., Phys. Rev. oi( liyjjp; 

^Millikan, R. A. and Neher, H. V., Phys. Rev. 47,204 (1935); 

_ c 50, 15 (1936). , . , 

16 Hess, V. P. and Graziadei, H. T., Terr. Magn. 4l, 9 (1936). 
17Forbush, S. E., Terr Magn. 42 , 1 (1937)* 

18 Schonland, B. F. J., Delatizky, B., Gaskell, J., Terr. Magn. 42 , 
137 ( 1937 ) • 

^Thompson, J. L., Phys. Rev. 54, 93 (1938). 

2 ?EpBtein, P. S., Phys. Rev. 537 862 (1938). 

21 Vallarta, M. and Godart, 0.. Rev. Mod. Phys. LI, 180 (1939)* 
2 fjanossy, L., Zeit. Phys. lo4 , 430 (1937)* 

28 Mallmfors, Arkiv. Mat. Astron. Fys. 30A, 12 (1944). 

2 ^MaHmfors, Arkiv. Mat. Astron. Fys. $2 , 8 (1945)* 

SJowight, K. Phys. Rev. 78, 40 (1950). 

26Alfven, H. Tellus 6, 23^-253, (1954) « 

2 ‘Dorman, L. I., Cosmic Ray Variations, State Publ. House for Tech. 
_ Q and Theor. Lit., Moscow, 1957* 

2 °Schluter, J., E. Naturforsch 62, 613 (1951). 

2 9piror, J., Phys. Rev. 94, 1017 (1954). 

30jory, F. S., P hys. Rev. 103 , 1068 (1956). 

3lLust, R., Phys. Rev. 105,~T827 (1957). 
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particles impacting vertically upon the earth, in order to 
relate 2 ^ 32,33 the local time (LT) of diurnal maximum (peak time) 
of CR to the assymptotic directions for the CR orbits which are 
responsible. From the CRDV peak times at a small number of ob- 
serving stations the direction of the axis of anisotropy was thus 
found 2 ^ to be a few hours before 18 hours LT and later found to 
be33 16.8 hours LT or^ 259° = 17*3 hours LT. Local time, LT, is 

also the longitude with respect to the earth-anti-sun line. 

35 

A study of the CRDV at Climax, Colorado and Huancayo, 

Peru showed that the daily variation (DV) in atmospheric temperature 
was not its principal cause, and this was later corroborated.-^ 

CRDV data from 1937 to 1951 were fitted to a differ- 
ential spectrum aE ^ times that of the omnidirectional CR for 
energies E^7-5 Bev spectral cutoff and a = 0 for E^7-5 Bev. 2 ? 

A review article on Cr 3^ discusses the CRDV. 

The daily variation bi-hourly counting rates (DVCR) for 
CR meson telescopes at Ahmedabad, India, showed^? that the peak 
times clustered in two groups at 03 hours and 11 hours LT at 
Solar minimum in 1954. These peak times Increased to 7-8 hours 
and 15 hours LT in 1956. Days when only the morning peak 
appeared occurred during a long-term decrease of daily mean 
intensity, and afternoon peaks were similarly associated with 

3 2 Nagashima, K.. Patnis, V. R. and Pomerantz, M. A., Nuovo Cimento 
„ 19, 292 (1961). 

^Duggal, S. P., Nagashima, K. and Pomerantz, M. A., J. Geophys. 

, Res. 66, 1970, (1961). „ , , 

. M., Phil. Mag. Vol. 6, #64, 573, (1961). 

Fonger, W. H. and Simpson, J. A., Phys. Rev. 94, 

, Adv. in Electronics and Electron Phys. 9, 129 (1957). 
37sarabhaiV,and Bhavsar, P. D., Supp. Nuovo Cim. 8, 299 (1958). 


-^Thompson, D 
35Firor, J.. 

1031 (1954) 
jDRose, D. C. 



increases. Sandstrom and Lindgren^® found that rejection of data 
for days of PD does not greatly affect the long term average 

2iO 

and phase of the reported CRDV. Kane J;7, has further discussed 

the affects on the CRDV of changes in the isotropic CR. Parsons 2 * 3 * 2 * 2 

found that even after FD are excluded, short term irregular CR 

variations depending on UT may alter the monthly average CRDV by 

not more than 0.11#. He attributed lack of simple agreement 

between peak LT's to some unremoved variation with UT. 

Directional meson telescopes at Uppsala, Sweden and 

43-45 

Murchison Bay, Norway, showed that n^ decreases with assymptotic 

latitude _A_, and is almost zero for -A. = 82°. Elliot's group has 

46-48 

done work on the CRDV. . The first and second harmonic CRDV 

amplitudes and peak times for all IGY neutron monitor stations 

have been calculated as averaged over the eighteen IGY months and 

4q 

plotted against both geomagnetic and geographic latitude. 

Although no curves are drawn through the points, a definite 
latitude dependence is displayed despite much scatter of points 
for stations of different longitudes. Messer schmldt-^ has pre- 
sented some of this Information. The CRDV has also been studied 
at low latitudes in India . ^ 


3®Sandstrom, A. E. and Lindgren, S., Ark. Fys. 16, No. 12, (1959). 
39Kane, R. P., Proc. Indian Acad. Sci. 52, 69-79 (i960). 
fOKane, R. P., Indian J. Phys. 35, 213 Tl96l). 

^Parsons, N. R., Tellus 12 , (4)7 I960. 

42parsons, N. R. J. Geophys. Res. 65, 3159 (i960). 

^3sandstrom, A. E., Dyring, E. and - Eindgren, S., Nature 187, 

1099 (I960). 

44sandstrom, Dyring, E. and Lingren, S., Tellus 12, 332, (i960). 
zSsandstrom, A. E. Am. J. Phys. 29, 187 (1961). 

4b Elliot, H. and Dolbear, D. W. N„ J. Atmos, and Terres. Phys. 1 , 
.,205 (1951). 

^'Elliot, H., Progress in Cosmic Ray Physics, 1 , p. 453* North 
.Holland Publishing Co., (1952). 

48EHiot, H. Phil. Mag. 5, 601-619 (i960). 

^9Schwachheim, G., J. Geophys. Res. 65, 3149 (i960). 
50Messerschmidt, W. , , Nat ur for s chung l3a, 734 (i960). 

5lRao, A. and Sarabhai, V., Proc. Roy. Soc. 263 , 101, 118, 127 
(1961). 
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The dependence of the CRDV upon longitude and month, as 
well as latitude and cutoff rigidity was presented by Cotten52 a t 
the April 1961 meeting of the American Geophysical Union, and was 
related to the locations of the several impact zones at each 
rigidity for near-equatorial assymptotic orbits, and to the 
higher multipoles of the geomagnetic field. Isoplot maps were shown. 

McCracken 5 ^ found m = 0 to 1 in a spectrum. 

Dattner and Venkatesen present some experimental results for the 
CRDV. 54 

Pomerantz, et al 5 f'^lve presented a theory for the CRDV, 
taking the amplitude as 

oO 

n l = -^2- /^S z (P,x)Jz(P)/ F ° j (cosA) n cos(^-^ E )dP 

\ '' P ' (eq.l 

where O^ o is a normalization constant, n = 1, and I, z, S z , J z , 
x, P, and P are as defined in section IV B on pages 30 to 33 . 

® OQ 

The assymptotic directions for CR orbits'^ Impacting vertically 
are the latitude with respect to the axis of anisotropy lying 
in the equatorial plane, and the longitude difference (^-^ E ) 
between the axis of anisotropy and the assymptotic orbit. This 
integral maps a CR anisotropy which is a cosine function of 
assymptotic longitude into a cosine function of longitude without 
consideration of the focusing into Impact zones accomplished by 
the geomagnetic field. 55 It yields average deflections by 

52 Cotten, D., J. Geophys. Res. 66, 2522, (1961). 

53MeCracken, K. G., Doctoral Thesis, Univ. of Tasmania. (1958). 
5^Dattner, A. and Venkateson, D., Tellus 11, 116 (1959). 

55Kelsall, T., J. Geophys. Res. 66, 4o47(T£6l). 



Prom, a least squares fit of fourteen 


the geomagnetic field, 
pairs of observed CRDV first harmonic coefficients to the values 
predicted by this theory, it was found^ that o(^= 0 for P<^7 BV/c 
and m = 0.4 for P^> 7 BV/c, with the axis of anisotropy at 16.8 hours 

LT. 

Forbush and Venkatesan^ 6 studied the yearly mean CRDV 
for 1937 to 1959 for ionization chamber data for Fredericksburg, 
Maryland; Huancayo, Peru; and Christchurch, England; and found 
that the yearly mean CRDV varied with a 22 year period. They 
established a statistically real CRD/2V but regard it, at least 
at H ua ncayo, as probably resulting from a systematic error due to 
friction at the barograph pen. 

Rao, McCracken, Venkatesan, and Katzman ^? -60 in their 
study of the CRDV at 29 stations, performed a harmonic analysis 
of the uncorrected neutron monitor DVCR and of the DV in atmos- 
pheric pressure at each station. They then performed a pressure 
correction upon. the uncorrected CRDV amplitude and phase, by a 
vector addition of the pressure correction of -.72^/mb corres- 
ponding to the amplitude and phase of the pressure DV. The 
resulting corrected CRDV values agree closely with those obtained 
by harmonic analysis of the cosmic ray nucleon counting rates 
which had been pressure corrected individually each bl-hour of 
each day. These investigators ascribe the second harmonic 
CRD/2V entirely to the second harmonic of the atmospheric pressure 


56p 0 rbush, S. E. and Venkatesan, D., J. Geophys. Res. 65, 2213 

(i960) . 

^Katzman, J., Can. J. Phys. 37, 1207 (1959)* 

5oKatzman, J. and Venkatesan, D. , Can. J. Phys. 38, 1011 (I960). 
^Katzman, J., Can.'J. Phys. 39 * 1477 (1961). 

6o Rao,u.,McCracken, K. G. and Venkatesan, D., J. Geophys. Res. 67, 
3590 (1962); 68, 345 (1963.). 



DV. They find that after pressure correction there is a residual 
CRD/2V at all equatorial stations except Huancayo, and none else- 
where . 

D. M. Thompson3 4 has studied the CRDV for neutron monitors 
at Makerere College, East Africa; Hermanus, South Africa; and 
Hwstmonceux, England for 1958 and 1959- No CRDV harmonic he finds 
could be caused by Improper pressure correction due to barograph 
pen friction because these stations use mercury barometers. He 
found a significant CRD/2V after pressure correction. 

Sternal j^ ag use( j power spectral analysis to determine 
the periodicities present in the CR neutron monitor counting rates 
at several stations and finds semi-diurnal, diurnal, 27 day, and 
annual periods. His analysis considers amplitudes only, ignoring 
peak times. It complements the harmonic analysis of other inves- 
tigators in that it does not assume a fundamental period of one 
day* yet finds one • 

Dessler, Ahluwalia and Gottlieb 10 have proposed a theory 
to relate the interplanetary magnetic field direction and strength 
to a CR anisotropy. Dattner and Venketesan 62 discuss many other 
possible causes for a CR anisotropy. 


^Stern, D., J. Geophys. Res. 
62 Dattner, A. and Venkatesan, 


67, 2133 (1962). 

d 7, Tellus 11, 239 (1959). 



I. D. Stormer Equations 


The equations of motion of a charged particle in a mag- 
netic dipole field were written and solved by Carl Stormer 3 and 
his followers since 1903 . These equations are shown in section 
VC , page 7& , of this report if Q is set equal to zero there. 

They have solutions in closed form only for orbits in the dipole 
equatorial plane. Numerical integration is required elsewhere. 
Because the field has symmetry around the dipole axis, the canonical 
momentum 2T conjugate to the magnetic longitude is a constant of 
the motion. No component of the ordinary angular momentum C X jp 
is conserved in general in a dipole field since the Lorentz force 
is not radial except for special orbits which are concentric circles 
in the equatorial plane. The magnetic rigidity P of a particle is 
defined as its momentum per unit charge. The Stormer constant 2 HP 
leads to the identification of allowed and unallowed regions, and 
to a lower limit for P called the cutoff rigidity P C (A, ,[* ) such 
that particles approaching from infinity will not impact upon the 
earth's surface at a magnetic latitude A from a zenith angle ^ 
and azimuth ^ unless their P ^ P c (A, ^ )• Por vertical 

incidence ^ = 0° and 

P c (A,0°) = P C (A) = 14.9 cos 4 A BV/c. (e 

Lemaitre and Vallarta^ 4 have modified the Stormer allowed cones 
to remove those directions whose orbits would have gone inside the 
earth at some other place. 


^stormer, C., The Polar Aurora , 

Lemaitre, G. and Vallarta, M. 

°^Lemaitre, G. and Vallarta, M. 

^Lemaitre, g. and Vallarta, M. 


Oxford, (1955). 
S., Phys. Rev. 43, 
S., Phys. Rev. 4$, 
S., Phys. Rev. 


87 (1933). 
719 (1936 . 
49, (1936). 
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Numerical integrations to determine the points of impact 
upon the earth for particles of various rigidities along orbits whose 
assymptotes at infinite distance (assymptotic orbits) have various 
directions have been performed by Stormer and his pupils, by 
several other workers for ^ = 0 31 an( j Kelsall^ for 
^ ^ 90°. This paper makes use of the impact zones55 which result 
from a broad parallel beam of non-interacting charged particles 
whose assymptotic orbits are at angle with respect to the dipole 

axis. If the beam came directly from the sun OC would be 90° at 

the equinoxes, except that the rotation of the earth carries the 
dipole axis around on a cone of half-angle 12 . In Kelsall 1 s-^ 
figures 6 it can be seen that for energies of less than 15 Bev 
there are two impact zones in each hemisphere, an early one at 
about 03 hrs LT, and a late one at 09 hrs, symmetrical about the 

equator. For energies higher than 15 Bev there is one zone 

centered on the equator but with essentially two centers . Both 
early and late impact zones show impact times later for high 
latitudes than for low or equatorial latitudes. When O<.^90 0 the 
impact zones shift away from symmetry. Kelsall's results show 
that the beam is focussed down from a broad area at infinity to 
much smaller areas of impact on the earth's surface. Focussing 
factors C (herein called f) are tabulated there. ^5 

I. E. Cosmic Rays in a Non-Dipole Geomagnetic Field 

Consideration of which orbits can reach the earth from 
infinity in a field model for the earth which includes some terms 
beyond the dipole terms has been made by approximately correcting 
the Stormer cutoff rigidities at vertical incidence (eq. 2) for the 
effects of the higher multipole terms. Cutoff rigidities appropriate 
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C*7 

for the eccentric dipole field model have thus been calculated. 

Quenby and Webber? provide formulae for a more correct 

set of P values which takes into account not only the dipole field 
c 

but also the effect of the local magnetic field, weighted according 

to the average relative importance of the various multipoles as 

averaged over the earth. The values of rigidity cutoff at the 

various IGY neutron monitor stations were computed according to 

6Q 

Quenby and Webber ' s formulae by Cogger? It has been shown that 
the earth's ring current cannot produce drastic reductions of 
apparent cutoff rigidity for near-equatorially incident particles. 


6?Kodama, M., Rondo, I. and Wada, M., J. Sci. Res. Inst. (Japan), 

6®§5gger, i 9 l]] Atomic Energy of Canada Ltd.-1104, Chalk River 

Ontario, CRGP-965. „ , TT . llo i. ai 

69Akasofu, S. and Lin, W. C., Trans. Am. Geophys. Union 43, 4b 1, 

( 1962 ). 



II. Data Examined 

T P 

Standard local production neutron monitor * data is 
examined for this study of the CRDV because it requires correction 
only for atmospheric pressure^* ^ at the monitor. Meson detectors 
for example require an additional correction for the (unknown) 
temperature and pressure distribution aloft. A network of forty- 
nine standard neutron monitor stations distributed widely over the 
earth gathered continuous records of counting rates during the 
eighteen months of the International Geophysical Year (IGY), July 
1957 through December 1958. The data, consisting of the number 
of counts per two hour interval, was submitted to the Japanese 
I GY World Data Center (WDC ). 3 Several stations performed pressure 
corrections, using their two-hourly barometer readings, and their 
local value for the pressure coefficient. The remaining stations 
sent their barometer data to the WDC, where the correction was 
performed using a common pressure coefficient of 0.96$/mmHg for 
all those stations. In effect, the WDC then expressed the 2-hour ly 
counting rates in tenths of percent of each station's average 
counting rate, by taking for most stations 


1000 In ( bl - hourl y counting rate 

V long term average counting rate 


Since the departures from the average are usually only a few per- 
cent and ln(l+x) & x for x 1, this gives the same number as 


^■Simpson, J. A., Ponger, W. and Trieman, S. B., Phys. Rev. 90, 

934 ( 1953 ). 

2 Arthur, W., The Cosmic Ray Increase of I960, (Ph.D. thesis. New 
York University) 

3cosmic Ray Intensity During the I GY, National Committee for the 
I GY, Science Council of Japan, Tokyo, (i960). 


would be obtained if the 2-hour ly counting rates were actually 
expressed as a fraction of the average. The data from two stations 
were expressed in simple percent of the station's average, by the 
stations themselves. The WDC then added together the percentage 
counting rates from the same two hour interval of each of the 
days in a month for which a complete 24 hour record was available, 
and divided by the number of complete days included in the month. 
This gives twelve numbers which represent the daily variation part 
of the cosmic radiation (DVCR) for an average day of that month. 
The averaging smooths out or reduces the effective value of any 

4 

fluctuations which did not occur at about the same time each ay. 


^Chapman, S. and Bartels, J. . Geomagnetism , Oxford, (1940). 



III. Data Reduction Performed 


A. Harmonic Analysis 

When an entry in any station's original remarks suggested 
it would be prudent, a day of data was eliminated and the remaining 
data was re-averaged. The bi-hour ly values of the average DVCR 
for each station and month were then analyzed to obtain the co- 
efficients in 

6 5 

DVCR(t) - y- « a + 2 a m cos m mi + 2 b m sin m mi (ill 

m=l IT m=l b 

where for i = 1, ...12, t = 2i-l is the Greenwich or "universal" 

p 

time, UT, in hours. 1 For this a standard routine using discreet 
sums was performed on IBM 70^ and 7090 computers.^ The coefficients 
were also expressed in the form 

6 

DVCR(t) - a + 2 n m sin m (t-JZf m ). 

° m=l 

The phase angles fS m were converted to local time LT in hours by 
Mni = K + Longitude. 

This phase time is 6/m hours before the peak time. For recognition 
of stations and months showing a strongly first harmonic CRDV, 


2 Willers, F. A., Practical Analysis , Dover, p. 3^5 (19^7). 
3 Program written by E. Me hr, NYU Research Division. 



the relative amplitudes 


_____ j ill — c_ j • t • v 

n l 

were calculated. A very few station— months of data were rejected 
because their higher harmonics m = 3, ...6, were thereby found 
excessive. A more stringent rejection scheme was attempted but 
was abandoned because it eliminated meaningful first harmonic data* 
merely because the stations 1 second or third harmonics were charac- 
teristically high. The sine and cosine coefficients for a harmonic 
dial 6 in local time were computed by 

a mLT = "m 008 ^mLT 
b mLT - "m 3ln *>mLT . 

A sample of the IBM 704 output listing is shown as Table I. Much 
of this output was also punched on IBM cards to form the input to 

the next computer program. Averages of n^, n 2 * flu/p* ^lUT* '^2!UT and 

0 were made for a great many groupings of stations and months. 
^2UT 

Sine functions of the same period are correctly averaged by a 
vector sum of amplitudes, taking into account their phase angles, 
divided by the number of entries, M. That is, 

n m =KlT 2 + ^mLT 
^mLT = tan 1 ^nLT 

a mLT 

£ 

See page 18 for reference. 



where 
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X mLT " 2 \l,T k ' m " 1 ' 2 * 

k=l 

where X = a mLT , b mLT , n^, or ^ mLT . 


(7) 


These, and the standard deviations 


a 

x 


77 v2 


Z N (X k - X k ) 
k=l _ K 
N - 1 


(8) 


were calculated on an IBM 650.^ Simple arithmetic averages of n m 

and 0 were also obtained. Samples of the input and output listings 
/ m 

are shown as Table 3T. 

A large Forbush-type decrease (FD) occurring in a month 
would decrease the daily variation ordinate for the bi-hour in 
which it occurred, and somewhat for the following hours, even after 
averaging over one month. The apparent DVCR would then not be the 
true DVCR. ^ The CRDV involves a peak-to-peak amplitude, or change 
in ordinate, of 0.6 to 1.0$. A typical PD of 6# would make only a 
0.2 % change in ordinate and therefore have a small but noticeable 
effect on ^ and 0^. Days showing an PD were not removed from 
the record. A few such days were removed from Zugspitze's data, 
and the harmonic coefficients were not affected sufficiently to 
seriously affect the maps and graphs which will be described. It 
was therefore considered unnecessary to reject all the data for 
those days. A linear increase in counting rate, such as is 

^Program written by A. Lucas, NYU Research Division. 

5carmichael> H., private communication. 
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Table XI 


Sample IBM 650 Input -Out put 
Vector Averaging of Harmonic Coefficients 
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KEY TO OUTPUT FORMAT 


IDENT. OF 

LA3T INPUT CARD 

At 

1*1 1 «r 

#..T 

n. lt 

0... 

in. 
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N 

A|UT 

6, ur 

IA„1 

A, 

O- B 1>rr 

cr I l IL ,| 
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N 
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|n,| 

O- A, LT 
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cr | n, | 
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N 
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K/n.l 

cr A t LT 

<r b. lt 

cr | n./n,| 



THE LAST Z DIGITS INDlC AT C DECIMAL POINT POSITION . 

NEGATIVE NUMBERS ARE INDICATED BV N , 9, CR , OR * FOLLOWING THE DIGITS. 



approximated by the exponential return to normal that occurs during 
the several days following an FD would affect^ the harmonic coeffi- 
cients by adding a saw-tooth waveform to the DVCR so that the 
apparent DVCR includes this saw-toothed waveform. For this reason 
the recovery periods following the selected removed days of FD were 
separately analyzed, and so were the remaining days in those months. 
There was no significant difference between their first harmonic 
coefficients and those obtained for the entire month. 

The fact that the DVCR ordinates represent averages over 
the separate two-hour intervals, rather than instantaneous values, 
reduces the first harmonic coefficients by° 


sin ( 2 12 ) 


2 (t ) 


= . 95 ^ 


so that all CRDV amplitudes n x found herein must be multiplied by 
( .954) -1 to obtain the true amplitude. This provides an increase 
which is very small compared to random errors. It will be accounted 
for in the end. Similarly, the second harmonic amplitudes are to 
be multiplied by (.826) 1 . 


III. B. Dependences Found 

The CRDV amplitude n ;L and peak time were found to depend 
on geomagnetic latitude A, figure 1, and less strongly upon longi- 
tude and month. In figures 1. to and 10 the small circles repre- 
sent vector averages for one station over several months, the large 


6 Chapman and Bartels. Geomagnetism , Oxford (19^0). 
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circles represent vector averages over several months for two or more 
stations, and the crosses represent arithmetic averages over the same 
groupings. The arithmetic average amplitudes are necessarily larger 
than or equal to the vector averaged amplitudes. The corresponding 
circles and crosses plotted closer together indicate smaller spread 
in the phase angle of the averaged sine functions. This closeness 
helps indicate that the data grouping is appropriate. The straight 
envelope lines on the graphs are at one standard deviation o, 
equation (&) , above and below the arithmetic average. This spread 
includes not only random errors, but also systematic dependences on 
month, and for figure 1 on longitude as well. Dark smooth curves 

show the relationships presented. 

The curves become smoother and more clearly exhibit the 
main features of the relationships when they are plotted against 
geomagnetic cutoff rigidity P c . Here P c is regarded as a parameter 
indicating effective position in the actual field. P Q Increases from 
the poles to the equator. Figure 2 shows n^ and 0-^ averaged over 
0.2 BV/c Intervals of P c as calculated for the eccentric dipole 
field model. ^ This combines the data for northern and southern 
hemispheres. The Quenby -Webber rigidities 8 are more representative 
of the actual field of the earth. Their use shifts the P q values 
of many points, and changes somewhat the grouping of the stations (Tab. Ill) 
into rigidity intervals, thereby altering the average amplitudes 
and peak times shown and making smoother curves, as in figures 3 
and result, especially when northern and southern hemispheres are 
treated separately. The P Q and n^ values at the six relative 
extremums of n 1 which are found as one proceeds from the Arctic 


^Kodama, M. 
^Quenby, J. 


, Wada, M., and Kondo, I., J.Sci.Instr.51, 138, (1957 ) « 
J., and Webber, W. R., Phil . Mag ., 4, 907 1959. 



TABLE III 


NEUTRON MONITOR STATIONS 

which have been averaged together for plotting at the same 


Quenby-Webber cutoff rigidity P c , and those 

plotted separately 

Northern (and Equatorial) Stations (P c ) 

Plotted at 

Thule (0.00), Resolute Bay (0.00) 

0.00 BV/c 

Murchison Bay (0.06), Helss Is. (0.07) 
Churchill (O.ll) 

0.1 

College 

0.48 

Deep River ( 0 . 87 ), Sulphur Mt. (0.94), 
Ottowa (0.96) 

0.9 

Mt. Washington (1.16), Uppsala (1.17)* 
Yakutsk (1.19) 

1.17 

Chicago 

1.54 

Leeds 

1.71 

Lincoln (1.99)* London (2.16) 

2.08 

Herstmonceux (2.30), Gottingen ( 2 . 38 ) 

2.34 

Climax 

2.77 

Weissenau (3-22), Zugspitze (3.33) 

3.3 

Berkeley 

4.01 

Pic du Midi 

4.30 

Rome 

5.04 

Alma Ata 

5.47 

Mt . Norikura 

9.13 

Makapuu Point 

11.03 

Huancayo 

14.18 

(Kampala) 

14.51 

Ahmedabad 

14.58 

(Lae) 

14.89 

Kodaikanal 

17.56 



Southern (and Equatorial) Stations 


Maws on 

0.57 

Mt. Wellington (Hobart) 

1.71 

Invercargill (Awarua) 

1.81 

Wellington 

3.20 

Sydney 

4.03 

Hermanus 

4.94 

Ushuaia 

5.89 

Buenos Aires 

10.70 

Rio de Janeiro 

11.47 

Mina Aguilar 

12.45 

( Kampala ) 

14.51 

Lae 

14.89 

(Kodaikanal) 

17.56 


The locations of the above stations can be found in 

reference (H-3) . All other stations listed there have been 

plotted on the monthly CRDV contour maps but omitted from the plots 

against P since they have data for very few of the included months. 
' w c 





the equator are shown in Table IV, Also shown are a single maximum 
and relative minimum in the southern hemisphere, where there are not 
enough stations to show as much detail as is found in the north. The 
persistence of these values despite considerable changes in the months 
which are averaged together, and whatever set of P c values is used, 
indicate that the two main peaks and ^troughs are real. The reality 
of these features should be questioned because it is possible to take 
a constant n x of 0.3# staying just within the envelope limits of + o, 
almost all the way to the Arctic where n-^ definitely falls steeply 
to 0.1# or less. Many prior investigators of CRDV have used a flat 
spectrum of this type. The detailed main features are not obvious 
until after some smoothing of data by taking averages such as have 
been described. Also, no stations have been eliminated by a severe 
data rejection scheme, as was done by several prior 9, 10,11 investi- 
gators. In this study, careful attention is paid to such details as 
the double peak of amplitude, and much new knowledge of the CRDV is 
obtained thereby. 

The peak time is seen to have a minimum near the equator, 

and to increase almost monotonically with latitude, or decreasing 

P . Strong dependence of n,, np and ,0 ilt u P° n the stations local 

position in the actual geomagnetic field as well as its latitude 

in the main dipole field has already been seen from the smoothness 

of their relationship to the Quenby-Webber P c , as compared with 

their relationship to the eccentric dipole P c . This local field 

dependence is better exhibited by isoplot (contour) maps of n^, 

n and 0 , figures 6, 9 and 7, as functions of latitude and 

2 1LT 

^Venketesen, D., Rao, U. R., McCracken, K. G., J. Geophys. Res. 67. > 
3590, (1962). 

10 Katzman, J., Venketesen, D., Can. J. Phys. 3°, 1011, (I960). 

11 Duggal, S. P., Nagashima, K., and Pomerantz, M. A., J. Geophys. 

Res. 66, 1970, (1961). 
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longitude. These CRDV maps have been made for all individual 
months, and for averages over several months. They all corres- 
pond very closely to isoplot maps of the horizontal or vertical 
intensity of the geomagnetic field at the earth’s surface. On 
figure 6 the double peak in the northern hemisphere and a simple 
peak in the southern hemisphere can be traced around the world, as 
indicated by the heavy black lines. The troughs are similarly 
traced by the heavy dashed lines. The equatorial flat of n-L shows 
a dark equatorial ridge line across India, Africa and South America, 
but a dashed equatorial trough line across most of the Pacific. 

The monthly n 1 maps show clearly a northern peak line and some- 
times indicate that it is a double peak. They show generally only 
a single trough near the equator, such as figure shows in the 
Pacific, but some months show the double trough at India and the 
Indian Ocean with a small equatorial maximum in between, as shown 
on figures S' and 6 . Figure 7 and all the other monthly or average 

mans of 0 show a minimum near the equator, with no well-defined 
* 'ILT 

maxima near the poles. All these ridge and trough lines are due 
to the main dipole field. Their departures from curves of constant 
magnetic latitude indicate the effect of non-dipole terms of the 
geomagnetic field. This is also indicated by ridge and trough 
lines which run approximately north-south on the CRDV maps, 
crossing the northern and southern peak zones and sometimes 
crossing the equator. The geomagnetic dipole equator is shown 
as a light dashed curve for comparison with the CRDV equators. 

The only serious ambiguities in CRDV first harmonic 
peak time occur in the Arctic where the amplitude is negligibly 
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small, and for several months at Hermanus, Capetown, South Africa, 
where a magnetic anomaly is located. Here the peak time, about 
03 hours LT, is almost twelve hours away from the peak time of 
neighboring stations. In drawing contours of peak time this must 
be shown as very early, or as very late. The choice between very 
early and very late was made in such a way as to obtain best agree- 
ment of the resulting contours with the contours for those months 
when no ambiguity arose. 

Figure ^ shows ridge and valley lines of relative maxima 
and minima of the Quenby-Webber cutoff rigidities drawn in geo- 
graphic coordinates . The maxima and minima which are oriented 
roughly North-South have been obtained at constant geomagnetic 
(dipole) latitude in the table in Quenby and Webber's paper, 8 and 
drawn upon the map of the Quenby-Webber rigidity contours taken 

from the same table. The P values would not show any such 

c 

extremums if the field were a dipole field. These zones are very 
similar to the ones on the CRDV maps for the months of symmetrical 
impact and most other months. This Indicates that the relative 
maxima and minima crossing the impact zones from North to South 
on the CRDV maps are due to and indicative of the non-dipole 
features of the earth's field. 

III. C. Second Harmonic CRDV 

The second harmonic CRDV amplitudes n^ are plotted 

against Quenby-Webber P in figure 10 . These values are obtained 

c 

after correction for atmospheric pressure. Northern and southern 
data are separated and shown "back to back" so as to display n^ 



from the Arctic through the equatorial region to the Antarctic. 

A north-south anisotropy is thus shown. In general, n g in the 

south is about twice n^ in the north. The equatorial maximum of *3$ 

n appears south of what is normally regarded as the equator, that 

is, either the geographic equator, the centered dipole equator, 

the geomagnetic dip equator, or the Q-W rigidity equator. Huancayo, 

on the centered dipole equator, at P c = 14.15 BV/c and n^ of only 

.045$, is definitely north of the equatorial peak of second 

harmonic CRDV. This is shown by the isoplot map of n g , figure 9 . 

Comparison with figure 8 shows that a conspicuous saddle point of 

the n contours coincides the only saddle point of the P contours, 

2 c 

and that Huancayo is the only equatorial station near that saddle 
point. It is therefore not surprising that Huancayo 's second 
harmonic is far less than that of the other equatorial stations. 

Not only does a true second harmonic CRDV exist after pressure 
correction, but its contours and extremum lines match those of the 
Quenby-Webber rigidity more precisely than any first harmonic CRDV 
contours and extremum lines in the equatorial region, especially 
regarding the equatorial maxima east of Kodaikanal, and the 
equatorial saddle point east of Huancayo, and the north-south 
ridge and valley lines, respectively, extending therefrom. 

HI. D. Annual Variation of CRDV Amplitude and Phase 

The mean sidereal day is about four minutes shorter than 
the mean solar day. Within a day therefore, a harmonic analysis 
based on either day as period would give about the same amplitude 
and phase, whether the actual period of the CRDV was one solar day 
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or one sidereal day, or whether both periodicities were present. If 
both periodicities are present they could be added vectorially to 
get the resulting observed CRDV since their periods are almost iden- 
tical. A variation whose period is one sidereal day has a vector 
representation which rotates at one revolution per year on a solar 
time dial. If it is added vectorially to the representation of a 
variation whose period is one mean solar day the resulting locus 
is a circle. The radius of the circle is the amplitude of the 
sidereal diurnal periodicity and its displacement from the origin 
Is the amplitude of the mean solar diurnal periodicity. 

In figure 11 the monthly change of amplitude and phase 
are shown for one year for North America. The polar plot shows 
an annual circling of the vector over the first 12 months. This 
Is consistent with a constant vector in local time plus a rotating 
vector which changes its phase, with the correct sense, through 24 
hours each year, and Is therefore constant In sidereal time. Simi- 
lar analysis has been carried out for many regions of the earth, Table 3C. 
Most regions confirm this apparent cosmic ray stellar diurnal 
variation. None deny its existence. The data has been smoothed 
to eliminate some monthly fluctuations for some regions of the 
world whose data only poorly showed the stellar diurnal variation. 

The result is that the stellar diurnal variation is shown more 
clearly in those regions. Fourier analysis in two directions 
gives the Lissajous figure of best fit as an ellipse, the semi- 
major and semi-minor axes of which correspond to amplitudes of 
0.13 % and 0.076#. Thus 0.08^ is an estimate of the amplitude of 
the stellar diurnal variation in cosmic radiation. Superposed on 
this approximate circling is a monthly fluctuation in amplitude 
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and in phase due in part to monthly change of the angle between the 

dipole axis and the incident beam, as indicated on the right side 

ft 

of the figure and in figure 8 in the paper by Kelsall. The fluc- 
tuations in total proton flux over the whole earth for a spectrum 
of e”" 5 closely resemble the observed monthly fluctuations except 
for small phase shifts which can be caused by the beam not coming 
directly from the sun, or anti-sun direction, and some irregular 
fluctuations, which may be statistical but which are due in part 
to monthly shifting of the impact zones relative to individual 
stations. There are also superposed random fluctuations in ampli- 
tude and phase of the apparent CRDV due to such events as solar 
flares and Forbush-type decreases in the cosmic radiation. 

A similar Lissajous analysis was performed to find a 
universal time component, which would be caused by oscillation 
of the impact zones due to diurnal change in the angle between 
the incident beam and the rotating dipole axis. The effect is 
masked by a longitude dependence due to the eccentricity of the 
earth's dipole field, but its , 1 % amplitude is discussed on page 83 . 

III. E. Relationship to Other Variables 

By plotting n^ and J^LT altitude of station within 

geographic regions of similar position in the geomagnetic field, it 
was found that the CRDV is independent of station altitude, or 
depth in the atmosphere. 

The diurnal variation (DV) in atmospheric pressure is 
related to but is not the cause of the observed CRDV. There exists 
at the top of the atmosphere over several stations an almost purely 
sinusoidally (first harmonic) CRDV. The not purely sinusoidal DV 

12 ~ 

Kelsall, T., J. Geophys. Res. 66, 4047 (1961). 
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in pressure, or thickness of absorbing air, masks that CRDV so that 
it is not observed in the raw data recorded at ground level. For 
Zugspitze the DV of the nucleon counting rates uncorrected for 
pressure variation did not resemble a sine curve. Neither did the 
pressure DV. Yet after the effects of that pressure variation were 
removed from the data the CRDV emerged as a reasonably smooth sine 
curve. It follows that while atmospheric pressure variation does 
not cause the CRDV, it does hide the CRDV, unless it is properly 
corrected for. Improper pressure correction will cause errors in 
the reported CRDV. 

Since the CRDV as averaged over all stations has its peak 

13 

time at 1400 hrs LT, it has been suggested by some that the CRDV 
might be due primarily to the diurnal variation in atmospheric 
temperature, which peaks at about the same time. This cause is 
j^ulscl out by the fact that the CRDV at equatorial stations peaks 
well before noon. Also if the effect were due to meteorological 
causes, the isoplots of CRDV might be expected to show some relation 
to weather maps, which they do not. 


13 


Haymes, R. C., private communication. 



Table 1 


Groups of Neutron Monitor Stations Analyzed 
to Investigate CRDV in Sidereal Time 




Date on which CRDV in 

No. 


ST is in phase with 

of 

Stations 

CRDV in LT during 

St/d* * 


July 1957- June 1958 
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IV. Explanation and Analysis 
A. General Explanation 

The general features of the CRDV dependence on A or P c can 
be explained by the impact zones for a beam of non-interacting parti- 
cles in excess of the isotropic cosmic radiation at infinity. The 
assymtotic orbits of such a beam would have to be at nearly 90° to 
the earth’s rotation axis to produce a CRDV which is nearly symmet- 
rical with respect to the equator. Furthermore there are several 
reasons to expect that such a beam would be nearly parallel to the 
ecliptic.^ The symmetrical impact zones 2 for a beam at C*s. = 90 ° 
therefore are correct for some months and then give the simplest 
and most illuminating picture of the world-wide A and time dependence 
of the CRDV. If the rigidity spectrum of the anisotropy beam falls 
off steeply from some peak value just above a low rigidity cutoff, 
as the spectrum for the isotropic cosmic radiation is believed to 
do, ^ then there will be two CRDV amplitude peak zones as observed 
in the northern hemisphere and two in the southern, at the rigidity 
which lies just above the spectral cutoff. Examine for example 
the impact zones for 5 Bev on Kelsall's figure 6a. These center 
at about 45° and 57° geomagnetic latitude and are at an early and 
a late impact time, respectively. The impact times of these zones 
would be 03 and 09 hours LT if the anisotropy came from the sun, 
which it does not, and if the difference between geomagnetic and 
geographic longitude is for this purpose ignored. A station, in 
the course of its daily rotation, will pass through these two zones 
and thus experience two pulses of cosmic radiation in excess of the 
usual average counting rates as shown in figure 13. These pulses 


^Dattner, A. and Venketesan, D., Tellus 11, 239, 1959. 
^Kelsall, T., Jour . Geophys . Res . 66, 4035, 1961. 
^Quenby, J. J*, and Webber^ W” Phil » Mag ♦ 4, 657> 


1959. 
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wlll be due to particles from two different rigidity intervals, for 
example, around 5 BV/c in the early zone at A = 45° and around 10 BV/c 
or more in the late zone at A = 45°. The observed CRDV does not 
peak between 03 and 09 hours LT because the beam does not come from 
the sun. The A dependence of the observed CRDV peak time is due to 
the fact that the impact zones are at later times at middle and high 
latitudes than at low latitudes, and to weighting factors for the 
early and late zones. These factors, Y e and Yj^ will be derived in 
the next section. The second harmonic CRDV can be partly explained 
by the existence of two impact zones, early and late, in each 
hemisphere. The shape of the pulses experienced as a station 
passes through the Impact zones, together with the fact that the 
early and late zones are not spaced exactly 12 hours apart, will 
give rise to higher harmonics. The non-dipole terms in the geo- 
magnetic field account for the local details of the zones of high 
or low amplitude. 

The earth's rotation varies ©<, thus causing the position of 
the impact zones to depend on Universal Time (UT). Thus a detector 
kept at constant LT (by not rotating with the earth) would experience 
a small CRDV dependent on UT, as impact zones for different rigidities 
passed over the detector. 
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IV. B. A Theory for Obtaining the Rigidity Spectrum of the CRDV 
The curves of CRDV amplitude versus P show the integral 
spectrum of particle momentum per unit charge responsible for the 
CRDV. In differentiating these curves to find the differential 
spectrum, account must be taken of rigidity-dependent neutron monitor 
counting efficiencies, S(P,X), as is usual for the spectral analysis 
of the total cosmic radiation, and also of a variable focusing factor 
f(P), and other factors imposed by the interaction between a colli- 
mated beam and the main field of the earth. For purposes of spectral 
analysis the main field will be taken to be that of a uniformly mag- 
netized spherical earth, a centered-dlpole field. However, the 
most accurate description of the integral rigidity spectra Is ob- 
tained by using the most accurate available values for P at each 
observing station. For this purpose the Quenby-Webber P c will be 
used. Where it Is necessary to compare with the P dependence of 
impact zones in a centered dipole field, comparison will be made 
at corresponding cutoff rigidities. If the anisotropy were per- 
fectly collimated, an unlikely situation, then the observed CRDV 
at individual stations would have a double-pulse form. The almost 
perfect sine form of the CRDV at several stations, notably Zugspitze, 
indicates that the anisotropy is not well collimated. An anisotropy 
which at remote distances from the earth varies sinusoidally with 
direction in momentum space^ can be written 

(j(P) + j 0 ( p ) cos 9 : ) dfl dP 

giving the flux of particles per unit time per unit area within the 
rigidity interval dP around P, coming from within the solid angle 

^Nagashima, K., V. R. Potnis, and Pomerantz, M., Nuovo Cim . 19, 
292-330, 1961. 


(IV-1) 



d-fl around the direction G f measured from the axis of anisotropy, 
as shown in figure 14 . J(P) is the differential spectrum for the 
isotropic cosmic radiation (usually regarded as the total cosmic 
radiation), and J Q (P) is the differential spectrum responsible for 
the CRDV. This distribution will now be normalized to agree with 
the flux j(p) in an equivalent collimated beam. 

/j( P ) + j o (p)cos0^) d/ldP ©y' • dA=^j(P)cos0 J +J o (P)cos 2 0^dadPdA 


is the number of particles per unit time withinthe rigidity interval 
dP passing through an area dA oriented along the axis of anisotropy. 
The total flux per unit rigidity interval, from all directions is 

therefore 



j(P) = J C (P) 

After traversing the dipole field this anisotropic flux will arrive 
at the top of the atmosphere in two impact zones in the northern 
and two in the southern hemisphere. The DVCR flux within these 
zones will be 


IJ j(P)dP 


where V is a unit step function 
station is within an impact zone, 
the average counting rate of 


which Indicates where and when a 
This gives a DVCR in excess of 



2ir r? r° 


y(0,t,x)=2 

z 



Pc( 9 '1-r) 


at depth x in the atmosphere. Here p is an appropriate function 
showing the distribution of particles over direction of incidence 
at the location 0. All interactions with atmospheric particles, 
and the relationship between the fluxes j of cosmic ray particles 
of various atomic numbers z at the top of the atmosphere and the 

3 

response of the detector at the earth's surface, are specified by 

the factors S„(p,_2L_). The local time t is regarded as a differ- 
zv cos ^ 7 

ence of either geographic or geomagnetic longitude since these 

differences are nearly the same. Contributions from particles 

reaching the atmosphere not directly over the station are included 

by the integration over the hemisphere above the horizon. This 

5 6 

integral is analogous to the usual ’ form for the total omni- 
direction cosmic radiation counting rate 


2tt J 

r r d C 


l(0,x)=z 

z 


3 z( P ’l^j) J Z < P > dP 3in ^ d /“- 

U 6 A>p) 


R 7 

Rather than integrate the above it is customary * ' to approximately 
relate the vertical flux per unit solid angle 


I v (e,x)=z / s z (p,x)J z (p)dP 

Z Jp c (e,}=o)=P c (e) 


5simpson, J. A ., Ponger, W. and Trieman,S.B. , Phys .Rev. 90/93^, (195*0. 
^Arthur, W., Doctoral Thesis, New York University, 1962. 

7 Brown, R. R., Nuovo Clm. 16, 956, (1957). 
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5 

to the omnidirectional counting rate by the Gross transformation 


I - 



( 6 ) 


where L is the mean free path of the cosmic ray particles in the 

3 

atmosphere. The vertical flux can be written 


I v (3*x) 


S(P,x) J(P)dP 


(7) 


where s(P,x) is a gross specific yield function relating the proton 
part J(P) of the total cosmic ray flux to the counting rate of the 
detector. 3 It is this S(P,x)j(P) which is tabulated and graphed 
in the paper by Quenby and Webber. 3 if P c ( >/*) were independent 
of ^ ,yu. and S were independent of x, <y the Gross transformation 
would be simply 

I = | dXl I v = 2 ttI v 


the term L/(L+x) is a reduction by atmospheric absorption and does 

not include the allowed cone® effects of P variation. In similar 

c 

fashion, the difficulties of integrating over the unspecified 
distribution p will here be avoided by making a Gross transformation 
to get 




y(8,t,x)=_L_ s 


r 


x+L 


y (p,e,t/(t))f(p,^(t))s z (p,x)j z (p)dp 


(0)=P J0^=°) 


( 8 ) 


8 

Lemaitre, G. ., and Vallarta, M. S., Phys. Rev. 50, 49, (1936) 
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The double-pulsed form of the DVCR is expressed by the explicit 
dependence of upon t. Within the early impact zone fj = 1 when 

p (e,^(t))(p<P 1 (0), t o (e,^(t))^t4 1 (0,oc(t)). 


Within the late impact zone 


l - 1 


when 


P p (^(t))<P<Po(0^( t ))> P)V0^( t ))' t 2 (0)<t(t 3 (0) 


Otherwise, IJ = 0. These P limits for oL = 90° are shown in figure 
12 as taken from table 6 in Kelsall's paper^ for ^ 45°. Compar- 

ison of this figure with a similar graph (not shown) for o< = 100° 
reveals that P x does not change much as ©<. changes from 77*5° to 
102.5° during one rotation of the earth, nor does P Q , except for 
an oscillation at the equator that is second harmonic in t, but 

p P . and P make a large oscillation. Examination of Kelsall's 

2 3 4 

impact zone maps reveals that the limits t Q and t^ make a large 
oscillation but and t^ do not change much as c< varies during 
a day. In the early zone, t Q ^t(t^ 


y e (9,t,x)= L s 
e x+E z 


p i<®> 


f f e' P > 
P o=< S > 


«C(t))S 2 (P,x) . J 2 (P)dP 


( 9 ) 


In the late zone, tg^t^t^ 




y-r (0, t,x) — L 2 
u x+l z; 


f T (P,^(t))S (P,x)j (P)dP+ MP,*(t))S g M/J^{ 




( 10 ) 
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At all other times DVCR(9, t,x)=y(0,t,x)=O . 

P = P or P , whichever is greater. Both y and y T vary only 
oc o c* c u 

slightly with t during their time interval where it is brief, since 
then does not change much. The pulses can therefore be approx- 
imated as rectangular by dropping the o<. dependence. This approxi- 
mation will be close except for the late zone at low latitudes 
where the pulse duration ^ ^ — t^t^ l ar 8 e * However, since the 
integrand for the late zone is less than the integrand for the 

early zone at most P values within their range, the error is small. 
The limits t Q and t 1 also change during the time interval between 
them. The early pulse begins at local time t 0 = t o (0,°^(t Q ) ) 
where t^ is the corresponding UT, and similarly ends at t-^ = 
t )). The effect is that the pulse duration - t 1 - t Q 

should not be obtained from a single impact zone map for a single 
value of , but t and t should be separately obtained. Since 
is only approximately known anyhow this small change will be 
ignored and both T e and T L will be obtained from Kelsall's 
tables 2 for e/- = 90°. An approximate y(e,t,x) consisting of two 
brief rectangular pulses of height y 0 and y L , duration TT e and 
X L and center to center time separation A^ is Fourier analyzed 
to obtain first harmonic coefficients 


r 


r e /2 


tl'P+ %. /z 


a-, = 


1 - 


y e sin t dt + 1 
7 r 


y L sin t dt 




v-rjt 


= -Y l sin A^ 


(ID 
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r 


U/; 


b-, = 


1 

7T 


V COS t dt + 

°e 


r-it+ r 7 i ' 

1 | y T COS t dt 


TT 






yL 


= Y + Y t cos 
e L 

where Y e,L = ^e.L 3ln |^e,L (9) 

7T 

and t, T, and ^ are In hours or in degrees, where 1 hour = 15°. 
The first harmonic amplitude of the CRDV Is given by 



c = / Y 2 + 2Y Y t cos A 1 /' + Y. 
1 J e e L 


( 12 ) 

(13) 


(14) 


It will now be shown that an uncollimated anisotropy of 
the same total flux gives a sine form CRDV of essentially the same 
amplitude as c 1 just calculated for the double-pulse form CRDV. 
the anisotropy J (P) cos 0* previously normalized, will be written 

j (P) sin 0 cos 0 = j(P,0,0) (15) 

o 

in transformed coordinates in rigidity (momentum) space such that 
the polar axis is parallel to the earth’s dipole axis. The angles 
9 and 0 are then the geomagnetic assymptotic colatitude and longi- 
tude of the particle orbit. A simpler anisotropy 


j(P,J0) = j E (P) cos ^ 


(16) 
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vfjLtJi equatorial assymptotic latitude will be treated. It ignores 
some spread in the latitude of the predicted CRDV amplitude peak 
zones. To normalize this to agree with j(P)> the total flux per 
unit rigidity interval through an area dA oriented along the axis 
x of anisotropy is obtained as 

oLTT 


r 


J x (P)<iA = 


j(P)dA 


j(P,0)-dA dp 




J (P) cos fidfi dA 
E 


or 


j(p) = ^(p) 


(17) 


Using j(P) = J Q (P) , (Eq 1 . 2)'. this' becomes 


•) E (P) - (f V P) 


( 18 ) 


so that the flux at large distance from the earth is 


j(P, / 0)dPd / 0 = 4 j (P) cos ft dP dp 

3 

within the increment dp in direction angle. From a broad beam at 
each assymptotic direction p a portion is selected by the geo- 
magnetic field and focussed to a pair of impact zones in each 
hemisphere for each P of (pulse) width and after deflec- 
tions if Q and ■h- At any local time t in geomagnetic coordinates 
a station will be within the late impact zones for assymptotic 
directions fc - ~ ?l) and simultaneousl y wlthin 

the early impact zones for^JZf^ - ^'f^p^pg + ‘5^'e) an< * rece l ve 


(19) 



at the top of the atmosphere a contribution 
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f ( p )5 J oz( p ) 003 ^ d ^ dP 


( 20 ) 


Prom the early zones over the station at t the total contribution 
to the DVCR detected at the surface is 


L Z 
x+L z 


r p(e) 


± f 0 (P)S f ,(P,x)j oz (P) cos 0 d)0 dP 


e'* ' z 




( 21 ) 


which is 


r 


(cos 


L Z 
x+L 2 


P,(e) 


8 f e (P)S z (P,x)3 oz (P) sin iT e (P,0)dP 


?*(*) 


( 22 ) 


From the late zones over the station at t the total contribution is 


L s 
x+L z 




_p 3 (e) M( 9 ) + i^ L ( p,0) 


.(e) 


5. f L (P)S z (P,x)J oz (P)cos 0 d)2f dP 




^(0)+iT T (P,0) 




Pa(e) 


i f L (P)S z (P,x)J oz (P)cos JZ) d0 dP 


^(0)-lT L (P,0) 


(23) 


which is 
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The total DVCR (0,t,x) Is the sum of these two contributions. The 

CRDV differential spectrum j (P) for each particle type z will be 

oz 

written as a fraction k z (P) of the differential spectrum J Z (P) for 
the total isotropic cosmic radiation. 




- 40 - 


with the 2 now under the integral sign. These two components have 

different amplitudes and phases but the same period (l day) and can 

be added vectorially to obtain the resultant amplitude c^ As in 

Quenby and Webber^ a set of coefficients K is introduced such that 

z 


J (P) - K J-. (P) , J(P) s J-,(P) . 
z z ^ - L 


(31) 


It is now assumed that the CRDV differential spectrum J oz (P) is the 

same fraction K(P) of the omnidirectional cosmic ray differential 

spectrum J (P) for all atomic numbers z, that is k z (P) = k(P). 
z 

Then 

k z (P)J z (P) = k(P)K z (P)j(P) ( 32 ) 

and 

2 S (P,x)k (P)J (P) = 2 K S (P,x)k(P)j(P) 

Z Z Z *7 " z 


= S(P,x)k(P)j(P) 

where 

S(P,x) = 2 K z S z (P,x) 


(33) 

( 34 ) 


is the gross S defined by Quenby and Webber^ to relate the proton 
part J(P) of the primary cosmic ray flux to the total counting 
rate of the detector. Then 


P e,L = f sln \T e L (P,0)f e>L (P)S(P,x)j(P)k(P). 


(35) 


The 2 has not been omitted, it has been formally performed, 
z 


Now 


let 


Pi(e) 


Y s _L_ 

e L+x 


P e dP 


P (©) 
oc v 


( 36 ) 
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Y 

L 


= L 



(37) 


These integrals can be written more compactly by noting that outside 
of the above integrated ranges of P,0,the factor sin ^/^(P ,0) in the 
integrand is zero. With this understanding 


Y 


e 



ao 



o o 



(38) 


(39) 


The law of cosines gives 

Cj 2 = Y e 2 + Y l 2 -aY e Y L cos(tt-A^) 
or 

C 1 = 

and from the law of sines 
s in £ = Yl sin 

yV e 2 +2Y e Y L cosa^+y l 2 (4o : 

which gives the time lapse £ from the center of the early zone to 
the peak time. 


/ 


+ 2Y Y t cos 
e L 




+ Y t 


W 
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This is identical to the result for the collimated beam, except 

that here the 2 sin * in Y is under the integral sign where it 

2 

properly belongs. 

Observed CRDV amplitudes n are given in this paper as 
fractions of the total omnidirectional cosmic ray counting rate I : 


n^e) 


ci(e) 

i(0) 


(4i) 


The factors which depend on colatitude Q in this analysis based on 
a dipole field model, depend also on longitude in the actual field 
o.f the earth. 

P (0) is a monotonic function so that 
c 

n^e) = n 1 (p c (0)) (42) 

is single -valued. Then 


dn^e) dP c = dn-L 


dP c (0) d0 


de 


(43) 


and 


nj? c (e))i(? c (e))=^^ 



(44) 


By use of the Gross transformation, eq. 6, 



this can be written 



2rt v (P c (9))ni(P c (e))= |fA dP ] 2+2 °° 3 ^ P e d ^ P L dP+ f/ FLdP ] 2 1 

Differentiating this gives an integral equation 



3 

This form is a bit more complicated than the usual technique in 
which the derivative with respect to P c of the known total (integral) 
cosmic ray spectrum is set equal to the differential proton spectrum 
multiplied by the specific counting efficiency S(P c ,x). If the 
sin i factor in F is regarded as dependent on 9 then 







+ f f b sln| f e>L (P,e)f e#L (P)S(P,x)j(P)k(P)dP 

where the first term is zero for the late zone integral since 

, lTC (P ,9) = 0 for all 9. The 9 dependence of sin iT'is only 
sin tL ' d 

very approximately known over much of the range of the variables. 
Regarding sin yf'Ti?) as independent of 9 amounts to squaring up the 
impact zones for each energy into rectangular blocks in latitude 
and longitude or local time extent. A similar squaring up in 
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local time and pulse height has already been performed in order to 
yield the expression for the CRDV amplitude, and hence a concomitant 
complete squaring off of the impact zones is reasonable for calcu- 
lation purposes • With that approximation 



d / ® sin (P)f~(P)S(P,x)j(p)k(P)dP 

W J 3 2 e e 

\ (*> 


= | sin ^(PpVP^SfP,,*)^)^! fcfP 1 ) 

sin §7>oc) f e( P oc )S(P oc' x > J < P °c>^°£ k(P oc ) 
3 e d0 


(48) 


where P = P or P , whichever is greater and 
oc o c 





T (P)f L (P)S(P,x) J Q (P)dP + d 



35/3 fJ'WAtf'b*’ 

r/e) 


--I slniT (P 2 )f L (P )S(P 2 ,x)J(P g )^2 lc(P 2 ) 

5 L L d d0 

+ X(\ ) £ fc) sft,*) Jfr,) * 

A© 

- | sin P 7 L (P !t )f L (P 4 ^)J(P„)^. k(P 4 ) . 


If, over a particular range of 0, as for A \56°, some of 
the limits P^, P^ do not exist, then for example, P^ and dP^/d0=O 

so its term drops out. Substitution of either of these forms for 
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the 0 derivative of the integral can lead to an iteration scheme 
for solving for k(P) where p c )> p 0 > but taking ^independent of 

0 is preferable where P c <C p 0 . 

At 53° the atmospheric cutoff rigidity 9 of 3 BV/c 
equals P^(0) so that Y e = 0 and the integral for c 1 becomes 


2TTl v (0)n 1 (0) = 





(50) 


and £ = Ay , which are exact relationships. 
Then, taking X L as independent of 0 


2irl v dn x dP c 

^ H5- 


+ 27m dl dP 
1 v c 

dP d0 
c 


-| slni?i,(P 2 )f L (I , 2 )S(P 2 ,x)j(lP | ):^ 


| Bln |T L (P 3 )f L ( P 3 )S(P 3' x)j(P 3 )dP 3 k(P 


cTS" 


- | sin ^T L (P 4 ) f L (P 4 )S(P V x)J(P 4 ) ^ k(P 4 ) 


(51) 


The variables I (P ), S(P,x), j(P), and P at the locations of the 

v c c 3, to 

stations are available from the papers of Quenby and Webber. The 
variables T(P). A'/W, P <) (8),P 1 (9), t^e), P 3 (0), P„(e) and f(P) 
are taken from Kelsall's tables 5 and 6. 


^Montgomery, D.J.X., Cosmic Ray Physics, Princeton, 19^9> p. 351* 
10 Quenby, J. J., and Webber, W. R., Phil, Mag. 4, 90, (1959). 
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f ( p) = I f o.=90 O + l f «c=100° 

is used as representing the average value of f on a day when the 
impact zones are most nearly symmetrical (o<»90°). dP/d0 = -dP/dX 
is obtained from the Stormer expression (eq. 1-2) for P^*) as 


dP c 

de 


sinX cos^X 


BV/c 
degree . 


( 53 ) 


Approximate values are available for every variable except the 
distribution function k(P), which is sought. 

Since the known terms in this integral equation are 
approximate, an approximate equation will suffice at least for the 
zeroth and first approximation solutions at X^53°. This is obtained 
by assuming in the expression for n-^ that Y e = Y^. This is reason- 
able over the range 0 o <^X<^50° , and is satisfied exactly at X = 4l.. 
“TK** 

c 1 = (2 Y Jl + cos 
Now replace 2Y ^ Y e + Y 

Lt 

and 

°1 = (Y e + V 003 2 


(54) 


( 55 ) 


c 1 (0) £ 


L 

x+t 


cos 




( 56 ) 
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Even if the above equality of early and late pulses is not approxi- 
mately met, the error in ^ incurred by using this formula is less 
than 10J* if Y A = 2 or | and less than 30 ^ if '* 1 / Y e * 100 or 1 / 100 
This approximation leads also to the simplification 

sin£ /y, 2Y sin 2 ^^ 

^ ij 

Y + Y 
e L 


or very roughly 


£ 



Y +Y r 
e J, 



for £ near to but these are not used for calculation. 

The amplitude relation becomes approximately 


2ttI (P (6))n (P (0))=cos 

V C X c 




( 37 ) 


and differentiating gives 

gri v fl* + 2irn l V". _ gTr Vl £ L ^_ l f . p +tf F Jp 

cos dP c d0 cosl^ip d p de (coslA^) 2 (tdl** del L 

2 Kie) 

(58) 

which can be solved for the unknown function k(P) contained in P 
by a straightforward Iteration process. For this purpose a zeroth 
approximation to k(P) la needed at p) 6 BV/c and is chosen so as to 
satisfy a further simplified form of this equation which is approxi- 


mately valid at low geomagnetic latitude. 
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Over the region 10^P c <(l8 BV/c, ^ is found to be nearly 
constant so the dn-j/dP term is zero. For P c <^12 BV/c, is nearly 

constant, so the d cos l/i^/d 9 term is zero also. Thus the two 
terms 


2ttI v (0) 

COS 1 (® ) 

2 ' 




(59) 


for 10<^ P c ( 9 ) 12 BV/c. This subsidiary equation has the solution 
n^(0)/cos = constant so that it is not necessary to 

separately assume n x and 4 are constant In order to assert that 


2im 1 (9) 
cos * A^( 0 ) 



(60) 


At these low latitudes P^ = P^ and it will be assumed for the 
zeroth approximation solution k Q (P) that P^ = P^ and P^ = P^. 

F and F_ can then be written as F(P) since the early and late 
e L 

zones don't have overlapping ranges of P. Then using e Q s . 57# 38* 39 
as 



equation 60 becomes 
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dl 


I v (0) dF 


or, substituting equations 35 and 7 


oO 




%g f 2 sin iT(P)f(P)S(P,x)j(P)k o (P)dP 


'2s in fSJlc 0 dP 

T 


p 0 (e) 


DO 


P (e) 
_£ 


■ OO 


d 


S(P,x)j(P)dP 


SJdP 


P 0 (S) 


p (e) 
c 


(61) 


which gives 

3ln |f(P 0 )f(P c )k o (P c ) 


3»i 


8 cos 

which is constant for all 10^P ^12 BV/c, so that there 


3 n i 


k (P) = 

° 8 cos 1 4l^f(P)sin 1 T 


>(62) 


k (P) ^ 0.004 


(63) 
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Pour forms divisible into two major classes are used to 

obtain k 1+1 (P) from k^(P) as solutions of the simplified equations 51&66, 

each valid over a separate interval of P c » The first class applies 

where Y = 0. It assumes Z is Independent of 9 or A. 

e 

At X)56°, P c <^1.4 &?/ c > 3^P 2 <V6 BV/c 


k (P ) = 
i+1 2 


3ir T dP dn + 3 ir n i SJ ( p c ,x 0 ) dP r 

“ ip l v c 1 1 c ° — < 


de dP„ 


d0 



+ sin \ T L (P 3 )f L (P 3 )Sj(P 3' X o ) ff3 k l (P 3 ) 

de 


V 


- sin | ^1<P 4 ) 

de 


J 


sin 


1 

2 


r L (P 2 )f L (P S )SJ(P 2' X 0>!!2 

de 


(64) 


is a solution to the exact equation, 51* 

For 56°^ ft and P 1 (e)^3 BV/c or k ± (P) = 0 for P^P 1 (e) j 

P >1.4 BV/c, P ^ 4.6 BV/c 
c / c. / 

then P^ = P^ and 

k(P)= - I v dP o f^l + ^ n l SJ(P o' x o ) 

de dP c 4 de 

sin \ T L ( p 2 > f L<y sj < p 2’ x o> 

de ( 6 5) 

is a solution to the exact eq./£l) which requires no iteration and 
is independent of the trial function k (P). These two solutions 
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are Joined by requiring continuity of the solution. The second 

class of solution is approximate and applies where Y g £ 0, at 

higher P . Again assuming "2" is independent of 9 or A, then 
c 

For 50.7°^ A, P c ) 2,i| 





- 77 T (P„)f T (P,)Sj(P„.x„) dp 

as 


S'L'VVa'^'V^' — ^ k.(p„) 

o n -L 




sin T e ( p 


)f (p )sj(p ) dP 

oc e oc oc oc 

de 


( 66 ) 


where for 50.7°} X )no°, 5.s)p,,)2- J t BV/c, 5. 2) P Q ) 3 BV/c 
then P oo = P 0 \P C , and for 40°)>X, P^ 5.2 BV/c then p o0 = p c )> P o' 


These two solutions are Joined by requiring continuity. Two 
alternate forms of the second class are valid if or ^ 
regarded as dependent on 9. These are ; for A ^40°: 


are 


k i = 



dP c d^ + if 1 m 1 sj(p #x) 


dP 


cos d0 dP c cosIa^ 1 


d0 

Equation (®7) 
continued on next page 
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- 3l I„n n < L - cos iA' 


TT v“l dA 


/ ( cos 

+ sln k (p )+ /"; 

d0 1 v4 /flV 




sin | T o (P c ^) f ^ P c )Sj(P c^ X o )dP c /d0 


e c 


(67) 


and 


k n ( p ) =, 

i+1 c 


-2ri v (P c (e)) dP c dn l + | T i 

— — — ao " - — ji *■ - 

COS t[ ^(e) de dP c cos ^&\p(e) 


( c^)5 


m 

+ / sin A 7< p . 0)f e (P)Sj(P,x)MP)dP 
^9 2 6 



+ 


/ 4.in£r<P.»)f L (PJ SJ(P)k.(P)dP 

■Ti(e) 


sin i 7” (P ,9)f (P )SJ(P ,x )dP< 

^ w e c ec oo 


d0 


These two forms are used to check the accuracy of the previous 
more approximate forms } equations 64, 65 and 66, numerically. 

The exact class of solutions k^(P 2 ) gives values over 
the range 3<^JP<^6.3, and the approximate class gives values over 
3 <^p<^n BV/c. Because of possible errors in the known functions 
yielding the "exact ^solutions (as well as errors in the functions 
yielding the approximate solutions), the approximate solution is 
not ignored over the range 3<^P<(6.3 BV/c for the purpose of 




( 68 ) 



- 53 - 

choosing the best fitting solution to the foregoing equations. 

That solution is used as the first trial for k^(P) in the iteration 
process for solving the exact equation for * 1+1 ( p oc ) for P c > ) 2,4 
BV/c, which is equation 46 in the form 



(Y el +Y L cosA^)sin iT 0 ( p ^) f A ( F oc )Sj(F oc* X o ) dP oc 

d0 

which applies when the F(P qc ) term^ P(Pg) term or P(P 1 )term. 


(69) 




Some choice must be made of months to be included in the 
spectral analysis. As will be discussed in Ejection Ta which con- 
siders the direction of the incident beam responsible for the CRDV, 
the impact zones are reasonably symmetric in the months September, 
October, November, March, April and May, when 90 • There 
exists a small but fairly consistent difference between Fall 

and Spring values of n_ and , as was discussed in paragraph 

1 HAT 0 

HID. For this reason six Fall months of o(~90 have been 
selected for spectral analysis, omitting the three Spring months. 

Only northern hemisphere data is used. Due to the scarcity of 
stations in the southern hemisphere and to the fact that few of 
them operated for all eighteen months, not as much detail shows ( ^ 4^ 
in the south as in the north. What detail does show there agrees 
well with the northern data. Not only was the CRDV data found to 
be independent of altitude, but it will be shown analytically that 
this independence is expected. For this reason data from stations 
of all altitudes enters the spectral analysis. The principal 
change which would be produced by the elimination of high altitude 
stations is further flattening of the broad equatorial region of 
almost constant n^. 
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The approximate equations (63) and then (66) which is a 
form of ( 58 ), followed by ( 69 ) which is a form of the exact equation 
(46) are used in solving the exact equation (45) in order to yield 
k(P) for P)>3 BV/c. These equations apply at P c ^ 2.4 BV/c, or over the 
main peak and the equatorial flat of n^. Equations (66) and ( 69 ) 
are valid iteration formulae provided that the term in P oc exceeds 
the terms in P g or P 1# Otherwise these terms are interchanged and 
the equation provides values of k(P) at Pg or Pi* not P Q or P c . 

Eq. 66 gives a solution which converges quickly, by the third 
iteration, and shows k(P)=0 for P^3*8 BV/c spectral cutoff, a peak 
from 3.8 to 5.5 BV/c, peaked at k( 5.*0 =0.0296, and low values 
ranging from k(6)=.001 to k(P^10)=.006 at p/6. This solution is 
inserted as the first trial k ± in the solution of eq. 69, which 
requires calculation of trial values Y e ^ and Y T ^ by eqs. 70* 71 • 

It is thus found that <j\ Y Ti for P c ^>6 BV/c, so that eqs. 50, 51 
apply there almost exactly, and give k(P) for 12. 9<^ ^ 17 • 7 Sv/'c. . 

They show k(P) is almost constant, k(P)~.004, at high P, with a 
small rise to .005 around P = 15 BV/c. Equations 66, 51 and 69 
indicate k( Picons tant at moderately high 10<P<20, and eq. 50 is 
then used to adjust that constant to the result k(P>6)=.0037 which 
gives the correct (observed) n^(P c ) at P c ^6 BV/c. Eq. 45 gives 
the same value. Equation 69 does not converge in an osculating 
manner, but oscillates. That is, if a trial solution k^ is put in 
which is too low, the resulting k 1+1 will be too high. In order to 
obtain convergence, therefore, the trial solution k^(P) is replaced 
by one intermediate in value between k^(P) and k^ ^(P), thus damping 
the oscillation. Six such Iterations produced convergence. The 
resulting k 1+1 (P) are integrated in eqs. 70*71 in order to satisfy 
eq. 45. In this way the solution converges to a peak from 3.8 to 
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5.7 BV/c, peaked at k(5.4)=.068, k(p)5.7)=.0037> and k(P)=0 for 
P<3.8 BV/c spectral cutoff, as shown in fig 16 . Equation 50 is 
a form of eq. 45 which applies at P c <(2.2 EV/c, over the small peak 
of n 1# and eq. 51 in the forms 64, 65 is used to solve it for k(P) 
for 3<P 2 <6.3 BV/c. It gives k(P)=0 for P^3*85 BV/c spectral cutoff, 
and shows a peak between 3-9 BV/c and 6.0 BV/c peaked at k(5«3)=0.0125, 
with very low k Just above 6 BV/c, as shown in fig 1 6 • Both the 
small peak and the main peak in n 1 (P (J ) thus result from a single 
peak in the differential spectrum k(P), located between about 3.8 
BV/c and about 6 BV/c, and peaked at 5-3 or 5.4 BV/c. However, the 
values from the small peak and from the main peak of n^ disagree by 
a factor of 5.38 as to the value of k at the peak. This discrepancy 
is largely due to uncertainty in the "known” ;f.unct ions A 'jj (0), 

i (p c ),T(p)> f ( p )> sj(p), p o ( 0) ' and p 2 (0) * The values of W' 

S j(p ), and dP^/dQ are especially poor for low P Q at the small peak. 
The discrepancy is also due to uncertainty in the values of n^(P c ), 
a for which is large, as shown in fig 5 . The height k discrepancy 
can be greatly reduced, to a factor of about 3* simply by taking 
the bend in n n (P„) at P = 5.4 BV/c to be less sharp. The peak in 
the solution k(P) of eq. 46 will then be less sharp, and lower. 

An attempt was made to take f e (P)^ p L (P) but to estimate them each 
from the values of f(P) tabulated by Kelsall together with the 
impact zone information given in his table 6, but the discrepany 
was then larger than when f e =f L =f was taken directly by interpolating 

(graphically) from Kelsall’ s table 5. 

For 7/P c <ll BV/c the values used for are guided by the 
almost constant small second harmonic relative amplitude, equation 
^9o) . At P c <^7 BV/c, can be obtained rather precisely from the 

sharply defined impact zones . Values of T are taken to have a 
simplified dependence on P, and no dependence on 0. 



The fact that one peak in k(P) leads to two peaks in 

n l(A(P c )) is due principally to the limits P Q , P^, Pg, P^ as functions 

of A, shown in figure 15 . At A <^38°, the k peak is outside the 

limits of integration of Y Q and When 39°< A <52°, corresponding 

to 5.5>>P ^ 2.2 BV/c, the k peak, and therefore a peak in F(P), 

eq. 35, is within the limits P and P of integration for the early 

o 1 

impact zone, thus producing large Y and the main peak of n . When 

e l 

52°^ A <'53°, corresponding to P 2*2.1, the k peak is outside the 

c 

limits of integration for Y g and Y L , in fact Y e =0, so ^ is small, 
producing the relative minimum n^ between the peaks. When 5^°/A <^62° 
corresponding to 1.9/> P ^ 0 BV/c, the k peak is within the limits Pg 
and P^ of integration for the late impact zone, thus producing 
especially large Y^ and the small peak of n^. Y^ is rather large 
at all A due to large ^ and tlie absence of an upper limit of inte- 
gration at A <^55°, hence the n^ peak due to large Y^ is smaller than 
the n x peak due to large Y e , which is zero or almost zero elsewhere. 

The k(P) peak is skewed, falling more steeply at its high 
P end, and so that the peak lies at the P value just below the P 

c 

value at the relative minimum n^ just south of the main peak. The 

spectral cutoff, k(P)=0, occurs at a P value P Q about 0.5 BV/c 

above the P value at the top of the main peak of n . Thus the 
C 1 

k(P) peak locations and spectral cutoff limits for periods other 
than September, October, and November, 1957 and 1958 can be obtained 
from table. nr. 

Although the integral spectrum n 1 (P (J ) can be perfectly 
fitted by a k(P) which is constant (k=00037) for p)> 6 BV/c, it is 
possible also to fit a k(P) which rises slightly from 6 BV/c to 
15 BV/c, and then falls off as 
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*P).A(a)\ P>15SE 
A = k(l5) 


(72) 


where m is very small. The constants A and m are chosen so as to 
yield the same values for Y e , and n^ as are obtained from the 
constant k. 

At high rigidity, using eq. 25* the CR flux (eq. l) can be 


written as 

(1+0.0037 cos e') J(P) dildP, P>6 BV 

c 


(73) 


or, using eqs. 25 and 72, as 


(l+15 m k(l5)P m cos 0 J ) J(P) dlL dP, P ) 1 3 BY 

(74) 


both in good agreement with the observed CRDV amplitudes. At lower 
rigidities the anisotropic part k(P) J(P) must be more elaborately 
specified, as shown in fig. , in order to describe the spike 

between 3.8 and 6 BV/c which is at least 3 and possibly 18 times as 

great as the anisotropy at p)> 6 BV/c. 

The small rise in k(P) near 15 BV/c was found after 
smoothing flat the fluctuations in n- L (P c ) near the equator, which 
are due partly to an unremoved longitude dependence, shown in fig. 
12. The k(P) found indicates o< = 0° in late December. When n 1 (P c ) 
is plotted for October 1957 through February 1958 or for December 
1957 through January 1958, then a pronounced rise in n-j/F c ) shows 
at 14,5 BV/c, suggesting a higher rise in k(P). Also the two peaks 
of n^ are then bifurcated or jagged. 
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The functions ^ sin S(P,x Q )j(P) F e> 

eq.(35) , used in obtaining the foregoing solutions are shown in 
figures 17a. and 17 t>, for low P. Integrations over them were per- 
formed graphically or by the trapezoidal rule. At high P> the 


integrations were performed analytically 3 using 

sin i T = 0.2 , P>6 BV/c (75) 

2 <- e 

sin |f L (P) = 0.03(P-5), 8<P 2 <25 BV/c (76) 

= 0.6 , P = ^ 2 S ? 25 BV//C (77) 

f(p) = 2.1-.0225P, 10 <'p<'49 BV/c (78) 

= 1 p>49 BV/c (79) 

S(P,x Q )j(P) = 200 P' 1 * 5 , P>10 BV/c (80) 


Because of the maxima in F e and F^ at 4.2. BV/c and 3*8 BV/c respec- 
tively, a simpler form such as 


k(P) = AP" m ( 8l ) 

with m = 0 , 0.4 , or 1 will also give two peaks in the Integral 
spectrum of hjLghtt and slope t different than those observed, but 
the spectral cutoff must be at P=3.85+«5 BV/e to agree with the 
location of the two peaks. This cutoff is lower than those found 
from the simpler theories. 
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IV. C. Altitude Independences of CRDV 

The principal dependence of both c^, (eq. 44), and I, 

(eqs. 6 and 5) upon altitude or upon atmospheric depth x occurs 
in the common factor (l+x/L), so that their ratio n^Ceq. 4l),is 
practically independent of x. Slight x dependence enters through 
the factor S(P,x) which Is integrated over different limits in eqs. 
5 and 44. Complete independence of x requires either the condition 

that 


( 82 ) 


with and P^ playing no role, so that the limits are the same in 
the two integrals, or else it requires the assumption that 


<fs(P,x) 

s(P,xJ 


1 

S(P,x) 


3s(p,x) 

7 x 


&x = K(x, 6x) 


which is nearly satisfied by Quenby and Webber's values for S . 
Then the change in vertical intensity 




Sl v = JL^ £x = /j(p) £ s(P,x)dP 

oo 

= /k(x, S x) j(P)S(P,x)dP 

J Pc 

= K(x, £x)l 


(83) 


( 84 ) 


Similarly, the change in n 1 is 
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Then, substituting the preceeding two expressions yields 


d x 


cos^ 

2irl 


K(x,Sx) Sx I" ^P e dP+^' L d p | + 


or 


2 ^- 0 . (86) 

*>X 

Thus the percent amplitude, 100 n^ of the CRDV, is expected to be 
approximately independent of altitude, insofar as the condition 
(eq. 83) upon S is met. The slight x dependence depends on the 
values of P-^ Pg, Pg and P^, which depend on 9. Hence the x 
dependence can be different for different groups of stations. It 
was found that n 1 increased slightly with altitude for some groups 
of stations and decreased slightly for others. The condition 
(eq. 82) for complete x independence is met in a centered dipole 
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field at X = 36°, P =6.5 BV/c. The almost complete independence of 

0 

the CRDV at all X permits spectral analysis using data from all altitudes, 
as was done, using X = for sea level in the function S(P,x) throughout 

the analysis. 


IV. D. An Explanation of the Second Harmonic of the CRDV 

A CR anisotropy which is purely first harmonic in assymptotlc 
longitude, such as eq. 19 > produces a DVCR, eq. 30, which is purely 

first harmonic in LT, despite the fact that there are two impact zones 
at each latitude. In order to have a second harmonic CRDV ( CRD/2V ) 
observed anywhere, the CR anisotropy must have an assymptotic longitude 
dependence which is more narrowly collimated than the cos jz5 dependence 
of eq. 19. Even then the CRD/2V will not be observed at all stations. 

A perfectly collimated anisotropy Includes all components 
in a Fourier analysis over assymptotlc longitude. This is the model 
discussed on pages 28 through 34, leading to two sharp rectangular 
pulses of DVCR, eqs. 9 and 10. For them the second harmonic Fourier 
coefficients are found, in a manner similar to that used for eqs. 

11-14, to be 


a 


2 




^4-T l A 


sin 2t dt 



dt 


= (l/V) sin T l sin 



and similarly 
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(l/ir) ^ sinT^ + sin T L cos (2 


Then, using = a 2 + b 2 and eq. 14 without the approximation of eq. 13, 
since is not small at low 



(y sin 7* ) 2 +2y e y L sin7^sin2^cos(2^)+(y L sin Tl) 2 

(y e sln | \ ) 2 +2y e y L sin| T e sin| T^os&f +(y^sin± T L )‘ 



(87) 


This gives n^/n^ = 1*1 a t X = 0° an ^ n 2^ n i = a® 

observed in the IGY data. Even when y & = 0, so that there is only 
one rectangular pulse, eq. 87 gives n^/n^ 5:5 ^ ^ or a P u ^ se infin- 
itesimal duration Tr • It is evident that n 
on the form of longitude dependence of the CR anisotropy, but also 
on the pulse shape and the pulse spacing A^of the associated CR 
pulses at the ground. Examination of the actual pulse shapes, as 
given by Kelsall's 2 table 6, indicates that at low X the pulse for 
50 Bev is broad and very nearly rectangular, and that the broad 
rectangular pulse approximation is not terrible at 25 Bev. A single 
broad rectangular pulse, of large 'Tr* does 8 ive n g^° low ** 
and a second pulse ((Y but at an appropriate & kJ/ gives an 
enhanced n^ The values of ^are again taken from Kelsall's 

table 6. At high X, ^«94° so that the second pulse makes n 2 much 
less than it would be for a single rectangular pulse. 




n^ depends not only 



At 39°<(A<(51 0 , within the main peak zone of n^, both 7^ 

and T t are sma ll so that e Q • 13 and 
L L 


(1/lr) y e,L ^"^e.L 53 Y e,L 


hold. At these A it was found in the solution of eq. 69 that 
Y «-Y t = Y, as in eq. 14, so that 


'J2 Y y 1 + cos ( 2 ) 


°2 * < y e + Y L> 


COS 




The relative amplitude than becomes 



The approximations in eqs. 14 and 88 were found to give errors less 
than 50$, and their ratio may be expected to differ from eq. 87 by 
even less. Eq.S9 does in fact agree very well with the observed 
n^/n^ at all A. It is to be remembered that eqs. 87 and 89 were 
derived for a perfectly collimated anisotropy, and must be modified 
for some other assymptotic longitude dependence. For example, 
n ^/n^ = 0 when that dependence is of form cos ^ . 

Another contribution to an enhanced CRD/2V at low A is 
suggested by eqs. 9 and 10. Diagrams such as figure \5 showing 
the P limits P , B , P_, P 0 , ...for oC= 90° and 100° clearly show 

C O 1 d 

a semidiurnal oscillation of P at low A, with P Q (0) unchanging 

at high A. It is not clear whether they also show a semidiurnal 

oscillation of P^, the lower limit for y , since the P limit 

2 
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diagrams are taken from Kelsall's table 6, which gives numbers only 

at 10 and 25 Bev, P at A~0 remains just below 25 BV/c. If it is 

2 

assumed that both the lower limits of integration P q and P^ have a 
semi-diurnal oscillation at very low A, this will produce a CRD/2V 
which has a fixed phase in UT and which will enhance the CRD/2V of 
fixed phase in LT which is introduced by the pulse shape and spacing. 

To show this take as found at low A from eq. 70 and the 

solution to eq. 69, and take 

P g = "F + b cos 2(t'-/3 ) (90) 

where t' = UT. Modification of eq. 10 gives 


y L = y + * P 2 


_d 

dP 


L f ( P , ) S (p,x )j(P)k(P)dP 
•z. / x+L o 


)) 


fi 


so that 


y = y + B cos 2 (t' - ) 


and 


_ 1 
a - — 

2 7T 


t/z 


y sin 2 t dt + — 

7T 


-r/z 


.r/z 

B cos 2 (t'-y# ) sin 2 t dt 


>-X/± 


Using LT 


( 91 ) 


(92) 


t = t ' + <p 
this becomes 

&2 = B sin(2 P+ 2/3 ) [Y- •i sin 47”^ /2tt 


( 93 ) 



Similarly 

b 2 = (y/ir) sinT + (B/2t r) cos(2^+2/?) [T + |f S in 4T)_ . 
Thus, if B = 0, then 

C 2 = (* Sln ^J // W 

which equals c for small *T but if B / 0 then c 2 can be larger or 
smaller, depending on the longitude at 7v^0°: 

c 2 2 = B S Uj S + I 3in 2 4t] + Tsm 4T cos (4 f+ 4/J ) j /4 tt 

+ (^F/ir 2 ) [T + | sin 4T]slnT cos (2f+ 2 /? ) + (y sinTA) 2 . 

The assymptotic longitude dependence of the CR anisotropy 
is not purely first harmonic, since a CRD/2V is observed at some 
equatorial stations. This is quite plausible, since cos gives 
a large angular spread for the beam. It is to be noticed that n^ 
is essentially the same for zero angular spread, eq. 14, or for a 
cos dependence, eq. 44, so k(P) as determined by eq. 45 or 69 
is unaffected by the conclusion that cos is not the complete 
assymptotic longitude dependence. 


IV. E. Application of LIouville ! s Theorem 

It Is commonly asserted. In studies of the isotropic CR 
and of the CRDV anisotropic CR, 10 that Liouville's Theorem 11 demands 


1 ^ > Ra6, U. R, , MfCraekehy K. G* and Venkatesan, ' D. , J. (Geopnys .Res . 

68, 345 (1963)." 

•^Goldstein, H., Classical Mechanics , Addison-Wesley, (1950). 
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that the CR flux at the top of the atmosphere equals zero for P<^P C , 
and e aual a the CR flux at infinity for P> Pc. 

If correct, this statement would forbid CR focusing by the geomagnetic 
field, and require f = 1 in eqs. 3 through 71. For a somewhat 
collimated anisotropic CR beam, the assertion Is not correct, and 

O 

focusing actually does occur as found by Kelsall. 

Liouville's theorem states that the density of points is 
constant in a 6N dimensional phase space where each point represents 
all the coordinates q* and momenta p of an entire system of N particles. 
The theorem can be applied to the isotropic CR case by taking each 
CR particle as an independent system, not interacting with Its 
fellows, so that the phase space Is six dimensional. Next the fact 
that p^ is constant for a charged particle in a magnetic field is 
Introduced in order to make an assertion about the volume in the 
three dimensional p s£>ace occupied by the N particles. The ordinary 
assertion for the isotropic CR case is that all directions of 
velocity are equally populated, so that for particles within dp 
around p the volume in p space is the constant 4 tt p 2 dp. 

Consider now an anisotropic CR particle beam at infinite 
distance from the geomagnetic field center, collimated so that not 
all directions of "p are equally populated. In fact, for a mono- 
directional, mono-energetic assymptotic beam such as treated by 
Kelsall, the volume in p space is zero, or if we treat an Interval 
dp around the constant "p, then it is an infinitesimal sphere 
4ir(dp)^/3. After deflection by the geomagnetic field, the particle 
momenta do not all have the same direction, in fact the beam is 
converging. The volume in p-space is then a much larger shell 
AilP 2 dp. Liouville's theorem requires that the volume in phase 
space be constant, so that the enlargement of the volume in p space 
requires a contraction in q space, thus increasing the particle 
density p In q space. Since v 2 is constant this increases the flux 



j*_^^ over the value it had at remote distances. Thus Liouville's 
theorem not only permits, but it demands geomagnetic focusing of an 
originally collimated charged particle beam. A similar discussion 
applies to the beam In a cathode ray tube, where magnetic focusing 
is also possible. Those papers on the CRDV which ignore magnetic 
focusing factors are making a substantial error at rigidities below 

6 BV/c. 
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V. Use of the CRDV as a Field Probe 
A. Assymptotic Direction of Incident Beam 

Information regarding the directions of the interplanetary 
(IP) and perhaps interstellar magnetic fields may be deduced from 
the direction of the axis of CR anisotropy responsible for the CRDV. 

This is the direction of the assymptotic orbit at the center of the 
beam, at a distance sufficiently remote from earth so as to be not 
yet significantly deflected by the earth's dipole field. To obtain 
this direction. A, the CRDV peak times are compared with the impact defec- 
tions \jj and and the predicted peak times 



theor 



6 + 12 . 


hrs 


corresponding to a perfectly collimated assymptotic beam from 12 
hours LT, perpendicular to the dipole axis (oC= 90° J. 1 

For the first approximation the above comparison was made 
using peak times from all months of the IGY, and without knowledge 
of the spectrum k(P). Eq. IV- 40 could then not be evaluated and 
£ sz—Aip was obtained by taking the centroid of the zone area at 

each energy in Kelsall's figure 6. The comparison 


(V-l) 


A = 12 + (p - 4 , 

X obs i theor 

then gave A » 20 hours, LT. The corresponding months of nearly 
symmetrical impact zones, when o<^90°, are November and May. 

For the second approximation, data might be taken from these three 
IGY months only. It is desirable, however, to average over more 


1 Kelsall, T„, J. Geophys. Res. 66, 4055* (1961). 



- 70 - 


than three months data to remove random monthly fluctuations. The 
zero-th approximation CRDV spectrum, k o J(F)=900O k 0 F" 2 ' 5 , has high 
values at low rigidity, so that Y e » Y L is possible, and then A « 24 
hrs LT. rnen equinoctial months would give the most symmetrical 
impact zones, and so they were included in the averages for the 
second approximation. For continuity of record the intervening three 
months October and April were also included. The ^ ersus 

tude, or versus P , curves obtained for these nine months of sym- 
metric impact agree closely with the curves for all eighteen IGY 
months. Therefore, the second approximation for A, using the same 
values of € as the first approximation, gave the same result. 

The differential spectrum k(P) stated on pages 56, 58 is 


based upon data from the months which would give almost symmetric 

impact, if the first approximate A is correct. That spectrum is 

used to evaluate eq. IV-40 and therefore the predicted peak time, 

eq. 1, at various latitudes, using the impact times \J/ from a zone 

map for oL = 9°° as drawn from Kelsall’s 1 table 6 for ^ 4 45 . 

This <£ theor is compared with the / 1LT curves for the same months 

using eq. 2 to obtain the third approximation for A. At most 

values of ? Q , € »A^and A = 17 to 17-6 hrs is obtained. However 

over the main peak of n_, Y g » Y L so 6 <kl^and the § theor are 

somewhat earlier than elsewhere. Eq. 2 applied to data in this 

region thus gives A > 18 hrs. If 6 is here computed 

using the Y from the k(P) peak of height as found from the 

e,L r- _ 

main peak region of n^(P c ), then the predicted <Q va ues 

over this region are too low (early) to agree with the observed 

^ 1LT* Better <^) agreement is obtained when a smaller peak 

height for k(P) is used, such as found from the small peak of 

n (P ). Thus the height discrepancy found in the peak of k(P), 
l v c 

page 56, is partially resolved if k(P) is required to satisfy 
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eq. IV-40 and V-l as well as eq. IV-45. 

North of 68° latitude not only is Y e = 0, but only the very 
high energy part of the late zone, P > P^, contributes. The peak 
time is then the not sharply defined centroid for high P ^ P^ and 
large , which occurs a little later than the low rigidity late 
zones. CRDV peak times in the Arctic are therefore exceptionally 

late. 

The fact that the source beam responsible for the CRDV may 

have a component from a direction opposite the sun's direction 

can be ascribed to a few related causes. They are the geomagnetic 

envelope due to the confinement of the geomagnetic 

field by the solar wind plasma, and the possibility that the sun 

casts a "shadow" for cosmic rays upon the earth. Our attention has 

focused upon the geomagnetic envelope. The significant feature is 

not the "tail" which stands away from the sun, but the "head" which 

faces the sun. 2 Here the horizontal component of the geomagnetic 

field is increased over its dipole value, while in the tail it is 

decreased. This increased field on the sunward side blocks out 

some of the otherwise isotropic galactic cosmic radiation, thus 

leaving an excess from the side away from the sun. 

That enhancement will cause a diurnal variation ctfi.thai-cut.bff 
rigidity at an observing station of form 


p (e)=p (e)<i+z r, (e)cos 

c c l i 1 


l(t-p ± )+Z R (Q)cos j(t*- 



(3) 


The UT=t 1 dependence is introduced by variation of as the earth 
rotates. Inserting this P n in eq.(lV-7)and using eq . (lV-80) 
gives 


^Suggested by A. Beiser, private communication 
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l v (e) = + 


= 4oo 


0 _p-°- 5 /0.5]p c 

T =’(el { 1+ l R k (9)0 ° S k (t ' 


- 1/2 


= 400 V(«) R 1 (S) ° 0S (t ' 5 R 2( e ) c032 ( t - 5 2 )+ -"} 

^ ( 4 ) 


This gives a diurnal and semi-diurnal variation in the observed 
cosmic ray counting rate I, due to no anisotropy in the cosmic 
radiation. If the ratio R^ is independent of Q or P c , then 
is independent of P Q . There is no reason to expect that R 1 (0) has 
a peak, and especially not a double peak, at high X. Thus the DV 
in P can not be entirely responsible for the CRDV, but it can make 
a substantial contribution. A CRDV vector due to DV in P c must then 
be added to the CRDV vector due to an anisotropic CR beam. This 
vector addition could produce a large shift In the apparent source 
direction A with only a small increase in amplitude n-^ and thus be 

responsible for the reported A ^ l8 hrs. 

A discrepancy between CRDV peak times observed at several 

P and peak times predicted from the k(P) presented can be ascribed 
c 

partly to the possibility that A depends on P. It is also partly 
due to an unremoved longitude dependence as shown in fig. 12, 
especially at Makapuu (P c = 11)* and Lae (P c = 14.9)* If A = 18 Hrs 
then oc « 90° on December 21. Little P c dependence is shown by 
CRDV peak times for October 1957 through February 1958 averaged, 
after the U.T. dependence is considered. 
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V. B. Earth's Magnetic Field 

The earth's quadrupole field is more significant than is 
indicated by the frequently stated 3 fact that its average value 
at the earth's surface is only 7$ of that of the dipole field. 

There are places on earth's surface near the dipole equator where 
the quadrupole field strength is nearly half the dipole field 
strength. At the dipole poles, the quadrupole field strength is 
almost zero, so its average strength is comparatively small. A 
similar statement applies to some of the higher multipole field 
strengths . 

On maps of the horizontal and vertical intensity (H and Z) 

of the earth's magnetic field, zones of relative maxima and minima 

of H and Z can be drawn. A CRDV relative maximum should occur 

where the dip is large, that is, wherever a minimum line for H 

nearly coincides with a maximum line for Z. There are two such 

places in the southern hemisphere and two in the northern, 

and they coincide the observed maxima of n^ for months of symmetric 

impact. This occurs because CR orbits, for rigidities not very 

much above P at high A, finish by spiraling tightly around geo- 
0 

magnetic field lines? Larger dip angles for field lines at a given 
A bring in particles of lower rigidity. Since S(P,x) J(P) peaks 
at low rigidity this provides larger I where the dip is larger. 
Since f (P)S(P,x Q ) j(P)k(P) is even more sharply peaked at low P, 
increased dip in the region of the n.^ peak zones will increase c^ 
faster than I, so that n = c^/l is largest where the dip is 
greatest. 

3 Quenby, J. J. and Webber, W. R., Phil. Mag. 4, 90, (1959). 



Similarly, A CRDV relative minimum should occur where the 
dip is small, that is, wherever a maximum line for H coincides a 
minimum line for Z. There are just two such places in the north 
and none in the south. These coincide the minima of n^ for months 
of symmetric impact. This shows that the CRDV amplitude is indic- 
ative of some of the details of the earth's magnetic field maps. 

Eqs. IV-14, 40, and 45 show that the CRDV peak time & , is 

1 1LT 

related to the amplitude n.^ through Y g and Y^. It follows that 
relative maxima and minima (ridges and valleys) on isoplot (contour) 
maps of ^ 1LT (Pig. 7 ) as well as of n (Pig. (o ) are indicative of 
similar relative extremums in the earth's magnetic field maps. 
Similarly, an isoplot map of the CRD/2V amplitude n^ shows ("Fuj \C>) an 
equatorial maximum which is Indicative of the geomagnetic equator. 

Most of the CRDV maps show four zones in the north. These 
are two North- South oriented ridge lines of relative maximum ampli- 
tude and also of phase, and two North-South valley lines of relative 
minimum amplitude and also of phase. Many of the maps also show 
four such lines in the south. These zones can be explained by a 
superposition of the earth's dipole plus quadrupole fields. The 
dipole field direction is from south to north near the dipole 
equator. The four quadrupole "poles" are situated near the dipole 
equator. These are the locations where the quadrupole terms of the 
scalar potential of the geomagnetic surface -measured field achieve 
their greatest magnitude, positive or negative. Above a positive 
such pole the quadrupole field is essentially radially up, thus 
making the field lines of the dipole plus quadrupole field dip 
more in the south and less in the north. This causes a CRDV 
relative maximum near a region south of the plus quadrupole pole and 
a relative minimum near a region north of the plus pole. Similarly 
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above a negative pole of the quadrupole the field lines of the com- 
bined field dip more in the north and less in the south. This causes 
a CRDV relative maximum near a region north of the negative quadrupole 
pole and a relative minimum near a region south of the negative 
quadrupole pole. Using this reasoning four quadrupole polar regions 
were obtained from the CRDV amplitude contour maps for November 1957 
and 1958 and May 1958. These regions will be called the "apparent 
quadrupole poles . " They agree within 15° with the locations of the 
quadrupole poles obtained from the spherical harmonic coefficients of 
the geomagnetic scalar potential. 

V. C. Orbits in Dipole Plus Quadrupole Field 
Cosmic ray orbits are governed by 
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The discussion of the preceeding section, TTB, ascribes the relative 

— * 

extremums of the CRDV to the dipole and quadrupole components of B, 
and so the n * 1,2 terms are most significant for describing the CR 
orbits responsible. The dipole term, in geomagnetic dipole coordi- 
nates, is 


V = a(a/r) 2 g° cos 9 . t 11 ) 

Because the quadrupole is so oriented that its poles lie near 9=ir/2 
in geomagnetic coordinates the dominant second degree term is m=2. 

This can be written 

V = a(a/rp h 2 sin 20 P« 2 (cos 0) ( 1 ^) 

2 2 

= 6a(a/r) 3 h g 2 sin 0 cos 0 sin 2 9 


by a suitable choice of 0 origin. In rectangular geomagnetic 
coordinates eqs. It, I 2 become 


V 1 = a3 S 1 ° Z// ^ x2 + y2 + z2 ) 3// 
V 2 = 6aV £ xy/(x 2 +y 2 +z 2 ) 3//2 


(13) 

(14) 


This dipole plus quadrupole magnetic field is then 


B x = ( 3a^g °xz - 6a^h^r )/ r^ + 

30a^h 2 x 2 y/r7 
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(15) 

B = (3a^g °yz-6a\ 2 x)/r 5 + 

y 1 2 

30a^h 2 xyz/r^ 

(16) 

B = (3aVy2-6aVx)/r 5 + 

2 -L c. 

30a^h 2 xyz/r^ 
2 

(17) 

where r = (x^+y 2 + z 2 ) 1 / 2 . 

Let 

dt = (ds)/v 

by making use of v = constant, eq. 

7, so that eq. 5 becomes 


P 2 

P d x/ds = B, z dy/ds - By dZ/ds 


(18) 

2 2 

P d y/ds = B dZ/ds - B dx/ds 
x z 


(19) 

P d 2 Z/ds 2 = By dx/ds - B x dy/ds 


(20) 


where P = mv/q 

It is convenient to combine the constants into the Stormer constant^ 
s k = y a 3 (-g°)q/mv' cm (21) 

and a coefficient 

Q' = 2ah|/g° jCm (22) 

when substituting eqs. 15, 16, 16 into eqs. 18, 19, 20. 

4 

Montgomery, D. J. X., Cosmic Ray Physics, p. 351, Princeton (1949). 
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These equations can be put in Stormer units by letting 
s=SS k 

X - xs k , y = YS k , Z « zs k 

r = RS k 
Q' = Q S k 


(23) 

(24) 


(25) 


Then, introducing 

1 = dX/dS 
m = dY/dS 
n = dZ/dS 

equations 17, 18, 19 become 


( 26 ) 

(27) 

(28) 


dl = 
dS 


-m[(2Z 2 -X 2 -Y 2 )/R 5 + 15 


QXYZ/R 


7 ] 


+n[3(YZ-Q^R 5 


+ 15 QXY 



(29) 


dm _ 
dS 


-n[3(xZ-QY)/R 5 + 15 QX 2 Y/R 7 ] 
+l[(2Z 2 -x 2 -Y 2 )/R 5 + 15 QXYZ/R 7 ] 


(30) 


dn = _i|3(YZ-QX)/R 5 + 15 QXY 2 /R 7 J 

+m[3(xZ-QY)/R 5 + 15 QX 2 Y/R 7 ] (31) 

2 . o 

with lengths and the coefficient Q in Stormer units. Taking h g /g 1 
_o.l then Q = 0.18 Stormer for P = 50.9 BV/c and Q = 0.06 Stormer 
for P = 5.9 BV/c CR particles. Equations 26-31 can be solved 



numerically for CR orbits * using 


df -a f(S+AS) - f(S) 
dS A S 

and successive approximations ^ (predictor-corrector method). 
Newton's method ^of extrapolation was found to give larger errors 
than successive approximations, where 

2 2 2 

err =1 + m + n -1. 

Some solutions indicate that orbits from zero assymptotic latitude, 
(z - 0, with n * 0,) are deflected by the quadrupole term so as to 
impact near the dip equator for this field model. Such orbits have 
been more extensively studied by Ruth Gall 6 using a similar set of 
equations . 

If Q is set equal to zero, the resulting equations for 
orbits in a dipole field are exactly those solved by Thomas Kelsall, 
using the Runge-Kutta technique and various initial conditions with 
n = 0 and R = 10 St'drmer. 


^Margenau, H. , and Murphy, G 
Chemistry, Van Nostrand, pp 
®Gall, R., J. Geophys. Res. 


. M., The Mathematics 

. 482-491, (1956). 

65, 3545, (I960). 


of Physics and 
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V. D. Effect of the Equatorial Ring Current 

The magnetic field above the ionosphere includes the fields 
of some systems of electric currents flowing in the neighborhood of 
the earth. A simplification of the field 7 of the equatorial ring 
current (ERC) will now be used to calculate the order of magnitude 
of the ERC 1 s contribution to the total deflection of CR particles 
by the geomagnetic field. In the region near the equator, the ERC 
field has opposite directions, approximately, inside and outside 
the current-ring. This decreases the field below the dipole field 
value over a region near the ring but inside it, r<r 2 , and increases 
the field outside the ring, r)r 2 » In this calculation it will be 
assumed that for r<r x the field is nearly the dipole value B 1# for 
r 1 (r<(r 2 it is B 1 b and for r^r 2 it Is B 1 a)> B,. This approxi- 

mation fits the values of Akasofu, et al 7 fairly well. 

The radius of curvature for a particle of rigidity P is 

R = P/B, so that deflection 



(34) 


3 

In the equatorial plane of a dipole field B = B c (r^/r) so the 
total deflection for equatorial CR orbits is 



2 * 

where (ds)^ = (rd0) + (dr) 


( 36 ) 


To get the order of magnitude of the produced by the 
ERC field let us take an orbit which is radially directed as it 
passes the ring, r g , and assume ds ^ dr for r^r^. Without the 

7 Akasofu, S.I., Cain, J . C . , and Chapman, S . , J . Geophys . Res .66, 4013,(1961). 
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ERC field the deflection suffered by a CR while traversing from r^ to 
r^ is then 



and with the ERC it is 



The deflections within the two regions in eq. 38 are almost equal 

and opposite. Various reasonable values for a, b, r , and r^, 

chosen to fit the approximate ERC field^ make (\^0]ring of eq. 38 

differ from Ad of eq. 37 by as much as 2° eastward or westward, for 

p = p =15 BV/c. This shift will be smaller for p)p c . Such a 
c 

shift will shift the direction A of the axis of CR anisotropy by the 
same few degrees, if unaccounted for in eq."2“l . Since A has been 
determined only correct to within about £0° the 2 shift can be 


ignored. 



Using a closer approximation to Akasofu, Cain and 
Chapman's approximate values for the ERC field. Dr. Untiedt 8 has 
more elaborately calculated the shift of the impact point of a 25 
Bev proton moving within the z = 0 plane with Stormer parameter 
_ _ 0.7, due to the introduction of the ERC field in addition 
to the dipole field. He found a shift of 0.8° towards the east. 

For z/Ohe has proposed a set of eight functions to approximate 
the ERC field within seven different regions of r. These are to 
be added to eqs. 15-17* 

It is probable that the shift caused by the ERC will be 
larger for CR of P <^15 BV/c which have non-equatorial impact points. 
For these the amount of the shift can be in latitude as well as 
longitude. It has been shown 9 that the ERC field does not drasti- 
cally alter the cutoff rigidity P c * Latitude shifting by the ERC 
is therefore small. Even at low P the longitude shifting can only 
show up in a shift of A. Any such error can be kept near to the 
plus or minus two degrees found at the equator by weighting the 
equatorial observed impact times more heavily than the high lat- 
itude ones found in the peak impact zone. This weighting will 
introduce errors larger than 2° however, since ^ g and ^ L are 
only roughly known at X = 0. 


8 Untiedt, J 
Gottingen, 
^Akasofu, S 


., Geophysikallsches Institut der Universitat, 
Germany, (private communication) 

., and Lin, W. C., Trans. A.G.U. 43, 461, (1962). 



V. E. Eccentricity of Geomagnetic Field 

The fact that the geomagnetic field is rather well rep- 
resented by an eccentric dipole manifests itself in that the mag- 
netic quantities, H, Z and field-dependent quantities P c , n^ / 1LT 
undergo one maximum and one minimum as one goes once around the 
equator. For n, , / nTm this is confused somewhat by the CRDV in 
UT, ^ which causes them to undergo a similar cycle as a function of 
longitude. 11 The CRDV at the equator was found to be representable (V 
by a vector of length n ^ .003 at constant / 1LT plus a vector of 
length £5 .0012 which circles in LT as longitude is varied, and in 
such a sense as to be fixed in UT. The vector of amplitude .0012 
is a vector sum of the CRDV vector fixed in UT plus an assymetry 
vector due to the eccentricity of the dipole. These rotate together 
in LT at a constant phase difference independent of longitude. The 
CRDV in UT may therefore be more or less than .12$. In the northern 
impact zone (at X = 50°) the circling is at times in such a sense 
as to deny the existence of a variation in UT, so the field eccen- 
tricity evidently dominates. 

V. F. Cosmic Ray Anisotropy Fixed in the Galaxy 

The regions of the world which show a small CRDV com- 
ponent vector rotating annually in LT in addition to a CRDV com- 
ponent vector fixed in LT are listed in Table IT. There agreement 
within one month is shown that the two vectors were in phase in 
late July during the I GY. African peak times indicate an earlier 
phase agreement, but this is to be disregarded since CRDV peak times 
at Hermanus, where there is a geomagnetic anomaly, often disagreed 
with expectations and with surrounding stations. 


^ 0 Kertz, W., Z. Geophysik 24, 210 (1959). 

^•^Kertz, W., private communication, Gottigen (I960) 
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If the vector rotating in LT and therefore fixed in 
sidereal time (ST) is a true indication of a CR anisotropy outside 
the solar system, then this anisotropic flux arrived at the earth's 
orbit from 14 hrs right ascension (RA) during the IGY. That is, it 
arrived from a direction in the galaxy which lay between the galactic 
center and the direction away from which the sun is moving due to 
galactic rotation. This direction is determined by taking A * 18 
hrs and the two flux vectors coming from the same direction on 
July 21, that is, from the direction to the sun on October 21, which 
is 13.7 hrs RA. A different value for A results in an equally 

different direction for the galactic CR anisotropy. 

It is possible that a boundary between the interplanetary 

and galactic magnetic fields exists with a "head" on the side toward 
which the sun is moving and a "tail" on the opposite side, similar 
in appearance to the geomagnetic envelope. 12 As discussed on 
pages 7/, 72, this could provide the excess CR flux from the tail 
side. This component would be added to a component perhaps from 
the galactic center or elsewhere 13 4o provide the anisotropy from 
l4 hrs RA. 

The eccentricity of the earth's orbit might conceivably 
produce the annual change in CRDV in LT, but it evidently does not. 

In July, when the two CRDV vectors are In phase so that the resultant 
amplitude n is largest, the earth is near aphelion. If the annual 
change in CRDV in LT were due to changing distance to the sun then 
n might be expected to be largest near perihelion. The perihelion 

and aphelion distances differ by only 4#. 

The obliquity of the earth's orbit with respect to the 

solar equator is another conceivable cause of the observed annual 


12 Beard. D. B. and Jenkins, E. B., J. Geophys. Res. 67 , 4895, (1962). 
13ross 1, B., Suppl . Nuovo Cimento 2, X, 275 (1955). 

■^Korff, S. A., Amer. Sci. 45, 293, (1957). 


change in CRDV. This effect might be expected to cause two maxima 
and two minima per year in n 1 . This would produce a small CRDV 
vector component which rotates semi-annually in LT, and not the 
annually rotating one. Such a semi-annually rotating vector of 
amplitude & .0005 is shown by several groups of stations. This 
very small effect may Indicate that the CR anisotropy in the solar 
system depends slightly on solar latitude, perhaps as a cosine 
function for example. 



VI . Cone lus I ons 


The amplitude and the local time of maximum (peak time) 
of the cosmic ray nucleon diurnal variation (CRDV) depend on latitude, 
longitude and month. Relative maxima and minima on their isoplot 
contour maps correspond to relative extremums in the geomagnetic 
field at or near the earth's surface, and in the geomagnetic cutoff 
rigidity. The CRDV map contours indicate that the geomagnetic 
dipole and quadrupole moments and the eccentricity of the geomagnetic 
field are significant in influencing the anisotropic part of the 
cosmic ray flux. 

When the amplitudes n ] of the CRDV of all northern and 
equatorial neutron monitors are averaged over groups of months 
within the I GY, and some pairs of stations with similar values of 
cutoff rigidity P c are averaged together to smooth the data, then 
two peaks of ni persist at approximately P = 1 and 4 BV/c. These 
are related to a single peak of k(P) in the differential spectrum 
(4/3) k(P) cos ft for the cosmic ray anisotropy. Here k(P) is a 
fraction of the isotropic cosmic ray flux at infinity and is 
the assymptotlc longitude with respect to an axis of anisotropy. 

P is magnetic rigidity, or momentum per unit charge. The k(P) 
peak lies between P = 3*8 and 5*7 BV/c, and k(P) = 0.0039 + 0.002 
for P^ 6 BV/c. At P^3.8 BV/c spectral cutoff, k(P) = 0. 

A meaningful semi-diurnal variation in cosmic radiation 
is found to be strongly dependent upon position on the earth. It 
indicates that the cosmic ray anisotropy is more narrowly colli- 
mated than a cos dependence. A significant part of the observed 
semi-diurnal variation is dependent on universal time. 

The CRDV is primarily caused by an anisotropy in the 
cosmic ray flux introduced by some mechanism within the solar 
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system. This anisotropic flux arrives from an assymptotic direction 
75 0 to 100° east of the sun, in the ecliptic plane. A small CRDV 
component dependent on local time is introduced by the increased 
geomagnetic field on the sunward side of the earth, so as to show 
an apparent excess flux from opposite the sun. If this component 
were subtracted, the direction of anisotropy stated above would be 
decreased. 

A component of the CRDV is dependent on universal time, 
with amplitude of about 0.1$. Another component of the CRDV of 
amplitude 0 . 1 $ appears to be dependent on sidereal time, and to be 
caused by a cosmic ray anisotropic flux in the galaxy which came 
from 14 hours right ascension during the IGY. 

It is possible to ascribe most, if not all, of the spatial 
variations in cosmic ray diurnal intensity variation on the earth's 
surface to the spectrum of the cosmic ray anisotropy and to the dipole 
and higher multipole geomagnetic field components. This suggests that 
insofar as the immediate vicinity of the earth is concerned, the ex- 
ternal geomagnetic field is fairly satisfactorily represented in its 
interaction with cosmic rays by the multipole coefficients obtained 

by analyzing the surface field. 

Several mechanisms appear responsible for the cosmic ray 
anisotropy. Further CRDV investigations should attempt to separate 
various sources and find their individual spectra and rigidity- 
dependent axial directions, rather than seek a single source 


direction. 
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CRDV AMPLITUDE nj AND PEAK TIME VERSUS 
GEOMAGNETIC CUTOFF RIGIDITY FOR THE 
ECCENTRIC DIPOLE GEOMAGNETIC FIELD 
MODEL 


24 r 


ZZ\ 


2Q\ 



THE LOCAL TIME OF DIURNAL MAXIMUM (CRDV 
PEAK TIME) OF THE COSMIC RAY NUCLEON COUNTING 
RATE PLOTTED AGAINST GEOMAGNETIC CUTOFF 
RIGIDITY P c as given by Quenby and Webber. The data 
shown, from northern and equatorial stations of all 
altitudes for September, October and November 1957 
and 1958, is typical of the entire IGY. 
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CRDV AMPLITUDES n 1 IN THE SOUTHERN HEMISPHERE 
PLOTTED AGAINST GEOMAGNETIC CUTOFF RIGIDITY P 
as given by Quenby and Webber. The data for each station 
is averaged over at least three of the months September, 
October, November, 1957 and 1958, and March, April, 
May, 1958. 
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averaged over September, October, November, 
1957 and 1958 and March, April and May, 1958. 
Contours are typical of the monthly maps of n^ 
for the entire IGY. 




Fig. 8. ISOPLOT MAP SHOWING RELATIVE EXTREMUMS 
IN THE GEOMAGNETIC CUTOFF RIGIDITY P c 
as calculated by Quenby and Webber, in geographic 
coordinate s 
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September, October, November, 1957 and 
1958 and March, April, and May, 1958. T1 
contours are typical of the monthly maps of 
ni for the entire IGY. 




Lissajous-Fourier analysis to determine the ampl 
and phase of a CR.DV fixed in sidereal time. Data 
shown is for 9 stations in North America from Jul 




Fig. 12. TYPICAL POLAR PLOT OF THE CRDV AMPLITUDE 
AND PHASE IN LOCAL TIME. Values shown are read 
from the amplitude and phase contour maps for September, 
1957, at uniform longitude intervals of 2 hours, along the 
CRDV equator. A Lissajous-Fourier analysis to find a 
CRDV fixed in Greenwich time is shown. 
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DOUBLE-PULSE FORM OF DVCR = y(t,0,oC) SHOWING 
QUANTITIES USED IN ITS HARMONIC ANALYSIS 
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Fig. 15. INTEGRATION LIMITS for rigidity integrals used in 
this analysis of the CRDV 
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Fig. 16. THE CRDV DIFFERENTIAE SPECTRUM k(P) 
given as the anisotropic cosmic ray flux in the 
direction of the axis of anisotropy, per cm^, per 
steradian, per second, per unit rigidity interval, 
divided by the isotropic cosmic ray flux. The 
anisotropic flux is taken to fall off as the cosine 
of the angle to the axis of anisotropy. 




Fig. 17. THE FUNCTION; 
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